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Accurate modelling of gravitational wave emission by extreme mass ratio inspirals is essential for 
their detection by the LISA mission. A leading perturbative approach involves the calculation of 
the self- force acting upon the smaller orbital body. In this work, we present the first application of 
the Poisson- Wiseman- Anderson method of 'matched expansions' to compute the self-force acting on 
a point particle moving in a curved spacetime. The method employs two expansions for the Green 
function which are respectively valid in the 'quasilocal' and 'distant past' regimes, and which may be 
matched together within the normal neighbourhood. We perform our calculation in a static region 
of the spherically symmetric Nariai spacetime (dS2 x S 2 ), in which scalar field perturbations are 
governed by a radial equation with a Poschl- Teller potential (frequently used as an approximation 
to the Schwarzschild radial potential) whose solutions are known in closed form. 

The key new ingredients in our study are: (i) very high order quasilocal expansions, and (ii) 
expansion of the 'distant past' Green function in quasinormal modes. In combination, these tools 
enable a detailed study of the properties of the scalar- field Green function. We demonstrate that the 
Green function is singular whenever x and x are connected by a null geodesic and apply asymptotic 
methods to determine the structure of the Green function near the null wavefront. We show that 
the singular part of the Green function undergoes a transition each time the null wavefront passes 
through a caustic point, following a repeating four- fold sequence S(a), 1/ttct, —5(a), —I/ttct, etc., 
where a is Synge's world function. 

The matched expansion method provides insight into the non-local properties of the self-force. We 
show that the self-force generated by the segment of worldline lying outside the normal neighbour- 
hood is not negligible. We apply the matched expansion method to compute the scalar self- force 
acting on a static particle on the Nariai spacetime, and validate against an alternative method, 
obtaining agreement to six decimal places. 

We conclude with a discussion of the implications for wave propagation and self-force calcula- 
tions. On black hole spacetimes, any expansion of the Green function in quasinormal modes must 
be augmented by a branch cut integral. Nevertheless we expect the Green function in Schwarzschild 
spacetime to inherit certain key features, such as a four- fold singular structure linked to the asymp- 
totic behaviour of quasinormal modes. In this way, the Nariai spacetime provides a fertile testing 
ground for developing insight into the non-local part of the self-force on black hole spacetimes. 



I. INTRODUCTION 



The last decade has seen a surge of interest in the nascent field of gravitational wave astronomy. Gravitational 
waves - propagating ripples in spacetime - are generated by some of the most violent processes in the known universe, 
such as supernovae, black hole mergers and galaxy collisions. These powerful processes are hidden from the view of 
'traditional' electromagnetic-wave telescopes behind shrouds of dust and radiation. On the other hand, gravitational 
waves are not strongly absorbed or scattered by intervening matter, and carry information about the dynamics at 
the heart of such processes. The prospects seem good for direct detection of gravitational waves in the near future. 
A number of ground-based detectors (such as LIGO pQ, VIRGO [2 and GEO600 [3J) are now in the data collection 
phase. 

Gravitational wave astronomy will enter a new era with the launch of the first space-based observatory: the Laser 
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Interferometer Space Antenna (LISA) [3]. It is hoped that this joint NASA/ESA mission, presently in the design and 
planning phase, will be launched within a decade. It will be preceded by a pathfinder mission, due for launch at the 
end of this year [5]. 

Black hole binary systems are a key target for gravitational wave (GW) observatories worldwide. Data analysis 
methods such as matched filtering may be applied to separate a weak GW signal from a noisy background [6]. An 
essential prerequisite for detection via matched filtering is accurate templates for the gravitational wave emission 
from black hole binaries. Breakthroughs in numerical relativity in the last five years have led to a rapid advance in 
the modelling of comparable-mass binaries, where the partners are of similar mass. Progress in numerical relativity 
continues apace. 

A key target for the LISA mission are the so-called Extreme Mass Ratio Inspirals (EMRIs): compact binaries in 
which one partner (mass M) is significantly more massive than the other (mass m). Mass ratios of \i = m/M > 10 -8 
are possible, for example for a solar-mass black hole orbiting a supermassive black hole [7 j. Mass ratios of up to 
m/M ~ 1/10 have been studied by numerical relativists |5]; smaller ratios are presently beyond the scope of numerical 
relativity due to the existence of two distinct and dissimilar length scales in the system. Perturbative approaches 
seem more likely to succeed in the extreme-mass regime. 

The smaller compact mass m distorts the curvature of the spacetime in which it is moving. Hence, rather than 
following a geodesic of the background spacetime generated by the larger mass M, the smaller mass follows a geodesic 
of the total spacetime [9]. However, if the mass ratio is extreme, the deviation of the smaller body's motion from the 
background geodesic will be (locally) small. The deviation may be interpreted as arising from a self-force, created 
by the smaller mass m interacting with its own gravitational field. To leading order, the self-force acceleration is 
proportional to m. With knowledge of the leading term in the self- force, one may model the evolution of the orbit 
and subsequent inspiral of the smaller mass, and compute the gravitational wave emission to high accuracy. However, 
finding the instantaneous self-force in a curved spacetime is not at all straightforward; it turns out to depend on the 
entire past history of the smaller mass, m. 

The idea of a self- force has a long history in physics. In the late 19th century it was well-known that a charge 
undergoing an acceleration in flat spacetime will generate electromagnetic radiation, and will feel a corresponding 
radiation reaction. The self-acceleration of a charged point particle in flat spacetime is given by the well-known 
Abraham-Lorentz-Dirac formula [10]. Radiation reaction implies that the 'classical' model of the atom (a point- 
particle electron orbiting a compact nucleus) is unstable. The observed stability of the atom remained a puzzle for 
many years, and provided a key motivation for the development of quantum mechanics. In the 1960s, DeWitt and 
Brehme [11] derived a formula for the self- force acting on an electrically-charged point particle in a curved background, 
and a correction was later provided by Hobbs [12]. The gravitational self- force acting on a point mass was found in 
1997 by two groups working concurrently and independently: Mino, Sasaki and Tanaka [13] and Quinn and Wald [13] . 
Shortly after, Quinn derived the self- force acting on a minimally-coupled scalar charge [15]. These developments are 
summarized in 2004/05 reviews by Poisson [16] and Detweiler [17 . In the subsequent period, a range of complementary 
approaches to the self- force problem have been developed p~8j [19] [20] [21] [22] . 

The self- force expressions for scalar, electromagnetic and gravitational cases take similar form [16]. In this paper, 
we restrict our attention to the simplest case: a point-like scalar charge q of mass m coupled to a massless scalar field 
$>(x) moving on a curved background geometry. The scalar field <&(x) satisfies the field equation 



where □ is the d'Alembertian on the curved background created by the larger mass M, R is the Ricci scalar, and £ 
is the curvature coupling constant. The charge density, p, of the point particle is 



where z(r) describes the worldline 7 of the particle with proper time r, g^ v is the background metric, g = det(g^), 
and S 4 (-) is the four-dimensional Dirac distribution. The field exerts a radiation reaction on the particle, creating a 
self-force [15] 



(□ - £R) $(x) = -4tt P (x) 



(1.1) 




(1.2) 



?self 



(1.3) 



which leads to the equations of motion for the scalar particle 



(1.4) 



where is the particle's four-velocity and is the radiative part of the field. Identifying the correct radiative field 
(which is regular at the particle's position) is the essential step in the derivation of the self- force [16] . Note that the 
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projection operator g^ v + u^u v has been applied here to ensure that u^f^ eU = 0. The mass m appearing in (1.4) is 



the 'dynamical' (and renormalized) particle's mass, which in the scalar case is not necessarily a constant of motion 
[T5] . Rather, it evolves according to 



dm 



(1.5) 



In other words, a spinless particle may radiate away its mass through the emission of monopolar waves. 

A leading method for computing the derivative of the radiative field, V v <&r, and hence the self- force, is based on 
mode sum regularization (MSR). The MSR approach was developed by Barack, Ori and collaborators [23| [24 ] [25] [26] 
and Detweiler and coworkers [T7J [27j [28] . The method has been applied to the Schwarzschild spacetime to compute, 
for example, the gravitational self-force for circular orbits [26] and the scalar self-force for eccentric orbits [29 . The 
application to Kerr is in progress [30j [31] . It was recently shown [32] that the gravitational self-force computed 
in the Lorenz gauge is in agreement with that found in the Regge- Wheeler gauge [T7J [27] . Further gauge- invariant 
comparisons, and comparison with the predictions of Post-Newtonian theory [33 , 34 are presently under consideration 

One drawback of the MSR method is that it gives relatively little geometric insight into the physical origin of 
the self-force. An alternative approach, based on matched expansions, was suggested by Poisson and Wiseman in 
1998 [36 . Their idea was to compute the self- force by matching together two independent expansions for the Green 
function, valid in 'quasilocal' and 'distant past' regimes. This suggestion was analysed by Anderson and Wiseman 
[37] , who concluded in 2005 that "this approach remains, in our opinion, in the category of 'promising but possessing 
some technical challenges'." The present paper represents the first practical attempt to implement this method. 

In the following sections we demonstrate that accurate self-force calculations via matched expansions are indeed 
feasible. We apply the method to compute the self-force for a scalar charge at fixed position on the product spacetime 
dS2 x S 2 (i.e. the product of a two-sphere and a two-dimensional de Sitter spacetime) introduced long ago by Nariai 
[38j[39]. We introduce a method for calculating the 'distant past' Green function using an expansion in quasinormal 
modes. The effect of caustics upon wave propagation is examined. This work is intended to lay a foundation for 
future studies of self- force in black hole spacetimes through matched expansions. The prospects for extending the 
calculation to the Schwarzschild spacetime appear good, although the work remains to be conducted. 

The remainder of this paper is organised as follows. In Sec. [II] we define the self- force and outline the Poisson- 
Anderson- Wiseman method of matched expansions. In Sec. |III| we consider wave propagation on the Schwarzschild 
spacetime. A radial equation of standard form is obtained via the well-known 'trick' of replacing the Schwarzschild 
potential with a so-called 'Poschl- Teller' potential. We show that a Poschl- Teller potential arises more naturally if we 
consider wave propagation on the 'Nariai' spacetime, whose properties are described in detail. 

Section IV is concerned with the scalar Green function on the Nariai spacetime. We begin in Sec. IV A| by expressing 
the Green function as a mode sum over angular modes and integral over frequency. We show in Sec. |IVB that 
performing the integral over frequency leaves a sum of residues: a so-called 'quasinormal mode sum', which may be 



matched onto a 'quasilocal' Green function, briefly described in Sec. IV C 



In Sec. [V] we consider the singular structure of the Green function. In Sec. V A we demonstrate that the Green 
function is singular on the null surface, even beyond the boundary of the normal neighbourhood and through caustics. 
We show that the singular behaviour arises from the large-/ asymptotics of the quasinormal mode sums. To investigate 
further, we employ two closely-related methods for converting sums into integrals, namely, the Watson transform and 
Poisson sum (Sec. VB). The form of the Green function close to the null cone is studied in detail in Sees. VD and 
VE[ and asymptotic expressions are derived. 



Section [VI] describes the calculation of the self- force for the sp ecific case of the static particle. For the Schwarzschild 

We show in Sec. 



VI B 



that the massive-field approach of Rosenthal 



spacetime, the static case has been well-studied 
may be adapted to the Nariai spacetime. This provides an independent check on the matched expansion calculation 
which is described in Sec. IVIC1 Relevant numerical methods are outlined in Sec. IVIDI 

In Sec. VII we present a selection of significant numerical results. We start in Sec. VII A| by examining the properties 
of the quasinormal mode Green function. In Sec. |VIIB] we test the asymptotic expressions describing the singularity 
In Sec. VII C we show that the 'quasilocal' and 'distant past' Green functions match in an appropriate 



structure. In Sec. 
regime. In Sec. VI 



D]we present results for the self-force on a static particle. 



We conclude in Sec. |VIII| with a discussion of the implications of this study. Throughout the paper, we employ 
geometrized units G = c = 1, and the metric sign convention { — h ++}• 
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FIG. 1: In the method of matched expansions, the tail integral is split into quasilocal (QL) and distant past (DP) parts. 



II. THE METHOD OF MATCHED EXPANSIONS 



Here we briefly outline the Poisson- Wiseman- Anderson method of 'matched expansions' [36, 37]. We start with an 
expression for the covariant derivative of the radiative scalar field [TSJ [16] , 

V„$* (z(t)) = -1(1 - 60lR^ + q(9^ + u^) Qd" + ^A« A ) + < ail (z(r)) (2.1) 

Here, R v \ is the Ricci tensor of the background metric and a v is the derivative with respect to proper time of the 
four- acceleration a v = The first two sets of terms are evaluated locally [14] [15] [16]. The final term <£^ ai1 is 

non-local; it is the so-called tail integral, 

< ail (*(r)) = q lim f ' V p G ret (^(r) ) z(r'))dr' (2.2) 

J— OO 

where Gretna?') is the retarded Green function, defined by 

B x G ret (x,x') = -An ] (2.3) 

together with appropriate causality conditions (which we describe in Sec IV A| ). 

Note that the tail integral depends on the entire past history of the particle's motion. Its evaluation is the main 
obstacle to progress. The tail integral (2.2) may be split into so-called quasilocal (QL) and distant past (DP) parts, 
as shown in Fig. [l] That is, 

c^(,(t)) = ^ L )(z(r)) + ^ DP )(z(r)) 

= q lim f T 6 S7^G re MT),z(T f ))dT f + q f T V ^G v Mr), z^dr' (2.4) 

e^0+J T _ Ar J_ OQ 

where r — At is the matching time, with At being a free parameter in the method (see Fig. [I]). 

The QL and DP parts may be evaluated separately using independent methods. In particular, if we choose At to 
be sufficiently small that z(r) and z(r — At) are within a convex normal neighbourhood [41], then the QL part may be 
evaluated by expressing the Green function in the Hadamard parametrix [42]. In other words, if z(r) and z(r — At) 
are connected by a unique timelike geodesic, then the QL integral is simply 

q- 1 ^ (z(t)) = - lim f 6 V,V(z(T),z(T f ))dT' (2.5) 

e^0+ J r _Ar 



wher e V(x , x') is the smooth symmetric biscalar describing the propagation of radiation within the light cone (see 
Sec. IV C for full details). The approach ultimately yields a series expansion for the QL self force in the coordinate 
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separation of the points x and x' . The Hadamard-expansion method is now well advanced for several spacetimes of 



physical relevance, such as Schwarzschild and Kerr [43] [44j |45j [46] H7J [48]. In Sec. IV C we apply this method to 
determine the quasilocal Green function and self- force in the Nariai spacetime. 

Evaluating the contribution to the Green function from the 'distant past' is a greater challenge, and is the main 
focus of this work. One possibility is to decompose the Green function into a sum over angular modes and an integral 
over frequency. In a spherically symmetric spacetime the Green function may be defined in terms of an integral 
transform and mode decomposition as follows, 

G ret (x,x f ) = — / dwe-™^ T(2l + l)PKcos 7 )^(r,rO (2.6) 



Here c is a positive constant, t and r are appropriate time and radial coordinates, and cos 7 = cos 6 cos 0' + 
sin sin 0' cos (<j) — (//), where 7 is the angle between the spacetime points x and x' . The radial Green function 
giuj( r , r ') ma y be constructed from two linearly-independent solutions of a radial equation. Since the DP Green func- 
tion does not need to be extended to coincidence (V — > r), the mode sum does not require regularization (though it 
may still be regularized if desired). However, Anderson and Wiseman [37] found the convergence of the mode sum to 
be poor, noting that going from 10 modes to 100 increased the accuracy by only a factor of three. 

In this paper we explore a new method for evaluating the L distant past 1 contribution, based on an expansion in so- 



called quasinormal modes. The integral over frequency in equation (2.6) may be evaluated by deforming the contour in 
the complex plane [49] [50] . This is shown in Fig. [6j In the Schwarzschild case there arise three distinct contributions 



to the Green function, from the three sections of the frequency integral in (2.6): 

1. A prompt response, arising from the integral along high-frequency arcs. 

2. A 'quasinormal mode sum', arising from the residues of poles in the lower half-plane of complex frequency uo. 

3. Power-law tail, arising from an integral along a branch cut. 

The three parts (1-3) are commonly supposed to dominate the scattered signal at early, intermediate and late times, 
respectively [49] [50]. (This may be slightly misleading, however; Leaver [49] notes that, in addition, the branch cut 
integral (part 3) "contributes heavily to the initial burst of radiation"). In this work, we investigate an alternative 
spacetime, introduced by Nariai in 1950 [38j [39], in which the power-law tail (part 3) is absent. We demonstrate 
that, on the Nariai spacetime, at suitably 'late times', the distant past Green function may be written as a sum over 



quasinormal modes (defined in Sec. IV B). We use the sum to compute the Green function, the radiative field and 
the self-force for a static particle. 

The key question addressed in this work is the following: how much of the self-force arises from the quasilocal 
region, and how much from the distant past? If the Green function falls off fast enough then only the QL integral 
would be needed, and, since the QL integral is restricted to the normal neighbourhood, only the Hadamard parametrix 
is required. Unfortunately this is not necessarily the case; Anderson and Wiseman [37] note that there are simple 



situations in which the DP integral in (2.6) gives the dominant contribution to the self-force. 

Usin g th e methods presented in this paper we are able to compute the retarded Green function and the integrand 
of Eq. \2.2\ as a function of time along the past worldline. We show that the DP contribution cannot be neglected. In 



particular, we find that the Green function and the integrand of Eq. (2.2) is singular whenever the two points z^{r) 
and z^{t') are connected by a null geodesic. We show that the singular form of the Green function changes every 
time a null geodesic passes through a caustic. On a spherically-symmetric spacetime, caustics occur at the antipodal 
points. 

On Schwarzschild spacetime, the presence of an unstable photon orbit at r = 3M implies that a null geodesic 
originating on a timelike worldline may later re-intersect the timelike worldline, by orbiting around the black hole. 
Hence the effect of caustics may be significant. For example, Fig. [2] shows orbiting null geodesies on the Schwarzschild 
spacetime which intersect timelike circular orbits of various radii. We believe that understanding the singular be- 



haviour of the integrand of Eq. (2.2) is a crucial step in understanding the origin of the non-local part of the self-force. 



As we shall see, the Nariai spacetime proves a fertile testing ground. 



III. SCHWARZSCHILD AND NARIAI SPACETIMES 



To evaluate the retarded Green function (2.6) we require solutions to the homogeneous scalar field equation on the 



appropriate curved background. In the absence of sources, the scalar field equation (1.1) is 



1 



A (y/^g^dv*) - £i» = 0. 



(3.1) 
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FIG. 2: Orbiting null geodesies on the Schwarzschild spacetime that intersect timelike circular orbits of various radii R = 6M, 
8M, 10M and 12M. The null geodesies are shown as coloured dotted lines, and the timelike circular orbit is shown as a black 
line. The spacetime point x is connected to xi,X2, etc. by null geodesies, as well as by the timelike circular geodesic. The 
Green function is singular when x — xi,X2, etc. Note that between R = 6M and R — 8M the ordering of the points X2 and 
X3 becomes reversed. 



For the Schwarzschild spacetime, the line element is 

ds 2 = -f s (r)dt 2 s + fsH^dr 2 + r 2 dSl\, d£l\ = dO 2 + sin 2 6d(j) 2 , (3.2) 

where /s(r) = 1 — 2M/r and the label 'S' denotes 'Schwarzschild'. Decomposing the field in the usual way, 

+oo +/] (S) ( \ 

u\J (r) 

^imus^imu s { x ) wnere <&i m(JJs {x) = " 

1=0 m=-l 



„oo +oo +Z u (S) / x 

= / du; s J2 E C im^im. s {x) where 9 lmu>3 (x) = ^J-Y lm {e, 0)e"^^, 



(3.3) 



where Y/ m (0,0) are the spherical harmonics, ci muJs are the coefficients in the mode decomposition, and the radial 
function (r) satisfies the radial equation 



r ^2 



+ «%-V l (S \r) 



«S(r)=0 



with an effective potential 



^ (S) W=/5(r) 



^ + 1) , /s(0 



2M N 



\ / Z(Z + 1) 



r / 



2M 



Here r* is a tortoise (Regge- Wheeler) coordinate, defined by 

fs r (r) r* = r + 2Mln(r/2M - 1) - (3M - 2Mln2). 



dr* 
dr 



(3.4) 



(3.5) 



(3.6) 



The outer region r G (2M, +00) of the Schwarzschild black hole is now covered by r* G (—00, +00). Note that we 
have chosen the integration constant for our convenience so that, in the high-/ limit, the peak of the potential barrier 
(at r = 3M) coincides with r* = 0. 
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Comparison of Effective Potentials for Schwarzschild and Nariai S 



Schwarzschild potential - 
Poschl-Teller potential ■ 



l = 10 




Radial Coordinate r. / M 



FIG. 3: Effective Potentials for Schwarzschil d \3.5\ and Poschl-Teller ^3. 
Schwarzschild tortoise coordinate is defined in (|3.6|) so that the peak is near r 

jo) 



radial wave equations. Note that here the 
= 0. The constants in (3.7) are Vb = 1(1 + 1), 



and a = 1/(V27M). 



A. Poschl-Teller Potential and Nariai Spacetime 



Unfortunately, to the best of our knowledge, closed- form solutions to (3.4) with potential (3.5) are not known. 



However, there is a closely-related potential for which exact solutions are available: the so-called Poschl-Teller poten- 
tial ED, 



(PT) 



(r*) = 



a 2 V 



(3.7) 



cosh (a(r* 

where a, Vo and r* ^ are constants (Vb may depend on I). Unlike the Schwarzschild potential, the Poschl-Teller 
potential is symmetric about r*°\ and decays exponentially in the limit r* — > oo. Yet, like the Schwarzschild potential 
it has single peak, and with appropriate choice of constants, the Poschl-Teller potential can be made to fit the 
Schwarzschild potential in the vicinity of this peak (see Fig. [3|. In the Schwarzschild spacetime, the peak of the 
potential barrier is associated with the unstable photon orbit at r = 3M. As mentioned in the previous sectio n (se e 
Fig. [2|, the photon orbit may lead to singularities in the 'distant past' Green function, and in the integrand of ( |2 2| ). 
Hence by building a toy model which includes an unstable null orbit, we hope to capture the essential features of the 
distant past Green function. Authors have found that the Poschl-Teller potential is a useful model for exploring (some 
of the) properties of the Schwarzschild solution, for example the quasinormal mode frequency spectrum [52j[53]. In 
this work, we hope to gain some insight into the 'Distant Past' integral on the Schwarzschild spacetime by using the 



exact wavefunctions for the Poschl-Teller potential, given later in Sec. IV A 1 

An obvious question follows: is there a spacetime on which the scalar field equation reduces to a radial wave 
equation with a Poschl-Teller potential? The answer turns out to be: yes [54, 55 j! The relevant spacetime was first 
introduced by Nariai in 1950 [5H1 [39]. 

To show the correspondence explicitly, let us define the line element 



ds 2 = -f(p)dt% 



f- 1 (p)d P 2 ^dn 



2> 



(3.5 



where f(p) = 1 — p and p G (— 1,+1). Line element (3.8) describes the central diamond of the Penrose diagram of 
the Nariai spacetime (Fig. |4|, which is described more fully in Sec. |IHB| Consider the wave equation (3.1) on this 



spacetime. We seek separable solutions of the form $(x) 
'Nariai'. The radial function satisfies the equation 



l ltJ N 



{p)Yi rn (6,(j))e- iuJNtN , where the label 'N' denotes 



V"" : 



du 



(AT) 



{u%-f(p)[l(l + l) + ZR])u£ ) N (p) = 



dp 



(3.9) 



where £ is the curvature coupling constant and R = 4 is the Ricci scalar. Now let us define a new tortoise coordinate 
in the usual way, 

dp* 



dp 



r\ P ) 



p* = tanh 1 p. 



(3.10) 
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Note that f(p) = sech 2 (p*) and the tortoise coordinate is in the range p* G (—00, +00). Hence radial equation (3.9) 
may be rewritten in Poschl- Teller form, 



dpi 



U 



cosh p* 



AN) 
1 Ioj n 



(P* 







(3.11) 



where Uq = 1(1 + 1) + 4£. We take the point of view that, as well as being of interest in its own right, the Nariai 
spacetime can provide insight into the propagation of waves on the Schwarzschild spacetime. The closest analogy 
between the two spacetimes is found by making the associations 



p* ^ ar*, tjsr^ats, cjn^^s/^ where a = l/(y/27M). 



(3.12) 



Fig.plshows the corresponding match between the potential barriers v£ S \r*) and V^ PT \p*) . In the following sections, 

we drop the label 'N', s o that t = t^ and u = uj>n. 

The solutions of Eq. (3.11) are presented in Sec. IV A 1 First, though, we consider the Nariai spacetime in more 
detail. 



B. Nariai spacetime 

The Nariai spacetime [38] [39] may be constructed from an embedding in a 6-dimensional Minkowski space 



ds 2 = -dZl + dZ< i 
of a 4-D surface determined by the two constraints, 



Zl + Zt + Zl = a\ 



Zi + Zt + Zi=a\ 



where a > 0, 



(3.13) 



(3.14) 



corresponding to a hyperboloid and a sphere, respectively. The entire manifold is covered by the coordinates 
{T, i/;, 6, (j)} defined via 



r 



Z\ = a cosh ( — ] cos^, 



r 



Zq = a sinh 

" a 

Zs = a sin 6 cos 0, Z4 = a sin sin 0, Z 5 = a cos 

with T G (—00, +00), i(j G [0, 27r), 6 G [0, 7r], <fi G [0, 2tt). The line-element is given by 



Z>2 = a cosh ( — ) sin'0, 



ds 2 



T 



-dT 2 + a 2 cosh 2 I ) dip* + a z dVt 



2 70 2 



(3.15) 
(3.16) 

(3.17) 



From this line-element one can see that the spacetime has the following features: (1) it has geometry dS2 x S 2 and 
topology R x S 1 x S 2 (the radius of the 1-sphere diminishes with time down to a value a at T = and then increases 
monotonically with time T, whereas the 2-spheres have constant radius a), (2) it is symmetric (ie, R^pa-r — 0), 
with Rp V = Ag^v, and constant Ricci scalar, R — 4A, where A = 1/a 2 is the value of the cosmological constant, (3) 
it is spherically symmetric (though not isotropic), homogeneous and locally (not globally) static, (4) its conformal 
structure can be obtained by noting the Kruskal-like coordinates defined via U = —(1 — AUV)(Zo + Zi)/2, V = 
— (1 — AUV)(Zq — Zi)/2, for which the line-element is then 



ds 1 



4dUdV 
' (l-AUVf 



dtti 



(3.18) 



Its two-dimensional conformal Penrose diagram is shown in Fig. [4] (see, e.g., [56]), where we have defined the conformal 
time £ = 2exp(T/a) G (0,7r). Its Penrose diagram differs from that of de Sitter spacetime in that here each point 
represents a 2-sphere of constant radius; note also that the corresponding angular coordinate ip in de Sitter spacetime 
has a different range, ip G [0,7r), as corresponds to its R x S 3 topology. Past and future timelike infinity V s1 coincide 
with past and future null infinity T ± , respectively, and they are all spacelike hyper surf aces. A consequence of the 
latter is the existence of 'past/future (cosmological) event horizons' [56, 571 EH] : not all events in the spacetime will 
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-oo-UV ■ 



X — i (t — — oo, p = oo; UV — j) X — i (t — — oc, p - 

FIG. 4: Penrose diagram for the Nariai spacetime in coordinates (ip, £). The hypersurfaces ip = and = 2tt are identified. 
Past/future timelike infinity coincides with past/future null infinity and they are all spacelike hypersurfaces. Thus, 

there exist observer-dependent past and future cosmological event-horizons, here marked as 7i% for an observer along ip = tt. 



be influentiable/observable by a geodesic observer; the boundary of the future/past of the worldline of the observer 
is its past/future (cosmological) event horizon. 

In this paper, we consider the static region of the Nariai spacetime which is covered by the coordinates {£, p, #, 0}, 
where p = atanh(p*/a) G (—a, -fa), p* = (v — u)/2 G (— oo, +oo), t = (v + u)/2 G (— oo, +oo) and the null coordinates 
{u,v} are given via U = ae _n / a , V = —ae v l a . This coordinate system, #,0}, covers the diamond-shaped region 
in the Penrose diagram (Fig. [4| around the hypersurface, say, ip = tt (because of homogeneity we could choose any 
other ij; = constant hypersurface). We denote by TL± the past cosmological event horizon at p = ±a of an observer 
moving along ip = tt; similarly, Ti\ will denote its future cosmological event horizon at p = ±a. Interestingly, Ginsparg 
and Perry [59] showed that this static region is obtained from the Schwarzschild-deSitter black hole spacetime as a 
particular limiting procedure in which the event and cosmological horizons coincide (see also [6Q | l6T | \62\ 163]). 

Note that there are three hypersurfaces p = 0, only two of which (those corresponding to ip — and 2tt) are 
identified (the one corresponding to ip = tt is not). Without loss o f ge nerality, we will take A = 1 = a. The 
line-element corresponding to this static coordinate system is given in (|3.8|). 



C. Geodesies on Nariai spacetime 



Let us now consider geodesies on the Nariai spacetime. Our chief motivation is to find the orbiting geodesies, the 
analogous rays to those shown in Fig. [2] for the Schwarzschild spacetime. We wish to find the coordinate times t — t' 
for which two angularly-separated points at the same c radius \ p, may be connected by a null geodesic. We expect the 
Green function to be singular at these times t — t' . 

We will assume that particle motion takes place within the central diamond of the Penrose diagram in Fig. [4j that 
is, the region — 1 < p < 1 (notwithstanding the fact that timelike geodesies may pass through the future horizons TL+ 
and Till in finite proper time). Without loss of generality, let us consider motion in the equatorial plane (6 = tt/2) 
described by the world line z fI (X) = [t(A), p(A), 7r/2, 0(A)] with tangent vector = [i,p, 0,0], where the overdot 
denotes differentiation with respect to an afBne parameter A. Symmetry implies two constants of motion, k = f(p)i 
and h = 0. The radial equation is p = ±i^(p 2 — Pq) 1 ^ 2 where po = y/l — k 2 /H 2 is the closest approach point and 
H 2 = h 2 + nw 2 . Here, w is the scaling of the afBne parameter and k = +1 for timelike geodesies, k = for null 
geodesies, and k = — 1 for spacelike geodesies. We still have the freedom to rescale the afBne parameter, A, by choosing 
a value for w. It is conventional to rescale so that A corresponds to proper time or distance, that is, set w = 1. Instead, 
we will rescale so that A = 0, that is, we set h = 1. 
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Let us consider a geodesic that starts at p = pi, cj) = which returns to 'radius' p — p\ after passing through an angle 
of Acj) (N.B. Acj) is unbounded, as opposed to 7 £ [0, n]). The geodesic distance in this case is s = —k(H 2 — l) 1//2 A(/>. 
It is straightforward to show that 



P{<t>) = Pi 



cosh (Hcj)- HAcj)/2) 
cosh (HA</)/2) 



hence 



The coordinate time it takes to go from p = pi , cj) 

Ati = 2/9*! + In 



= /9isech(iJA0/2). 

= to p = pi, cj) = is 

Vi-Po + V(l-Po)(Pi-^o) ^ 
,Pi + - \/(l-P§)(Pi-/»o), 



where p* 1 = tanh 1 (pi). Substituting (3.20) into (3.21) yields Ati as a function of the angle Acj), 



/ 1 - pisech 2 (i7A^/2) + tanh(tf A0/2W1 - p 2 sech 2 (#A(/>/2) \ 
Ah = 2p* x + In V 

\1 + Pl sech 2 (HAcj)/2) - tanh(#A(/V2)y 1 - p 2 sech 2 (#A(/>/2) / 
This takes a particularly simple form as p\ — > 1, 

Ati - 2^*, + In (sinh 2 (#A(/>/2)) , for pi -> 1. 
As A0 — > 00, the geodesic coordinate time increases linearly with the orbital angle Acj) 

Ati ~ 2p* x + HAcj), for Acj) -> 00, pi -> 1. 



(3.19) 



(3.20) 



(3.21) 



(3.22) 



(3.23) 



(3.24) 



In other words, for fixed spatial points near p = 1, the geodesic coordinate times Ati are very nearly periodic, with 



period 2ttH. Results fl3.19| ) , d3.22[ ), ( |3.23D and ( |3.24D will prove useful when we come to consider the singularities of 
the Green function in Sees. IV Dl and IV El 



IV. THE SCALAR GREEN FUNCTION 



Retarded Green function as a Mode Sum 



The retarded Green function for a scalar field on the Nariai spacetime is defined by Eq. (2.3), together with 



appropriate causality conditions. As described in Sec. [TT] the Green function may be defined through an integral 
transform and a mode sum, 



-1 r + OO + ic +OO 

G ret 7) = / ^^^(p*,p'J(2Z + l)P/(cos 7 )e- lw(t - t ' ) 

^ J — OO + ic j = q 



(4.1) 



where p* and t are the 'tortoise' and 'time' coordinates in the line element (3.8), c is a positive real constant, t — t' 



is the coordinate time difference, and 7 is the spatial angle separating the points. The remaining ingredient in this 
formulation is the one-dimensional (radial) Green function giu{p, p') which satisfies 



dpi 



U 



cosh 



(4.2) 



The radial Green function may be constructed from two linearly-independent solutions of the radial equation (3.11). 
To ensure a retarded Green function we apply causal boundary conditions: no flux may emerge from the past horizons 
HZ and (see Fig. Uk. To this end, we will employ a pair of solutions denoted uf^ and u^, in analogy with the 
>e. Thes< 



Schwarzschild case. 



rese solutions are defined in the next subsection. 
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1. Radial Solutions 



The homogeneous radial equation (3.11) may be rewritten as the Legendre differential equation 

d / x x " 2 

dp 

where 



(i - p 2 ) 



dp 



13(13 + 1) 



P 



1 



■p 1 







p = ±iu>, /3= -1/2 + 
A = ±v/(l + l/2) 2 + d, d = 4£-l/2. 



(4.3) 



(4.4) 
(4.5) 



We choose p = ioo, A = ^/(Z + 1/2) 2 + d and note that the choice of signs will not have a bearing on the result. The 
value of the constant £ in the conformally-coupled case in a D-dimensional spacetime is: (D — 2)/(4(D — 1)). Note 
that for conformal coupling in 4-D (£ = 1/6) the constant is d = 1/6, and for minimal coupling (£ = 0) we have 
d = —1/2. For the special value £ = 1/8 we have d = 0. The possible significance of the value £ = 1/8, the conformal 
coupling factor in three dim ensions, was recently noted in a study of the self- force on wormhole spacetimes [64] . 

The solutions of Eq. (| 4.3[ ) are Legendre functions of complex order, which are defined in terms of hypergeometric 
functions as follows (Ref - [65 Eq. (8.771)), 



i 



r(i- M ) \i-pj 



i 



■A 



/z/2 



2^1 



-P,/3 + 1;1 -m 



In the particular case p = 0, the solutions belong to the class of conical functions (Ref. 
the pair of linearly-independent solutions to be 



(in) 
/up) 



(p) 
(p) 



r(i -[i)p%(-p), 



r(i -fi)p^(p). 



(4.6) 

Eq. (8.84)). We define 

(4.7) 
(4.8) 



These solutions are labelled "in" and "up" because they obey analogous boundary conditions to the "ingoing at 
horizon" and "outgoing at infinity" solutions that are causally appropriate in the Schwarzschild case [50 j. It is 
straightforward to verify that the "in" and "up" solutions obey 



(up) +iuop* 



as 

as p* 



-00 
-00 



(4.9) 
(4.10) 



To find the asymptotes of near p = 1, we may employ the series expansion 



2 F 1 (a,b;c; z) 



T(c)T(a + b-c) 



T(a)T(b) 

T(c)T(c-a-b) 
r(c - a)r(c - b) 



(1-z) 



c—a—b 



E 

k=0 

&kbk 



(c-a) k (c-b) k (1 - z) k 



(c + 1 - a - b) k 
(l-z) k 



Lfc=o 



(l + a + 6-c) A 



where (z)^ = T(z-\-k)/T(z) is the Pochhammer symbol. In our case a 
It is straightforward to show that 



/in) 



(a 



4( out )> 



ifc! 

— 00, 

+oo, 



where 



A' 



j (in) 



(out) 



r(l - iw)T(-iuj) 
T{\+P-tw)T(-P-iuy 

r(i - iw)r(iw) 
r(i + /j)r(-/?)' 



fc! 

(4.11) 

l-/i and 1-z = (l-p)/2. 

(4.12) 

(4.13) 
(4.14) 




FIG. 5: Penrose diagrams for IN and UP radial solutions. 



with (3 as defined in (|4.4|). The "up" solution is found from the "in" solution via spatial inversion p — > —p; hence 



A (out) e 



— 00, 

+oo. 



The Wronskian W of the two linearly-independent solutions and u^(p*) can be easily obtained: 



dp* luJ dp* 
The one-dimensional Green function giui(p*,p'*) is then given by 

!,(p*K?(p'*)> P* < P'*> 



(4.15) 



(4.16) 



9iu(p*,p'*) 



w \ uZ(p*)u£( P :), P * > a, 

= \T{1 + (3 - p)T{-(3 - p)P%{- P< )P£{ P> ), (4.17) 
where p< = min(p, p') and p > = max(p, //). The four-dimensional retarded Green function can thus be written as 

^ oo 

G iet (x,x') = —^(21 + l)fl(co8 7 ) 



i=0 



/+oo+jc / 1 \ / 1 

dwe -4 "^-* } r - + i\ - iuj r - 



iA - ioj) P^ /2+iX (-p<)P^ /2+iX ( P >) (4.18) 



B. Distant Past Green Function: The Quasinormal Mode Sum 

As discussed in Sec. [ilj the integral over frequency in Eq. ( |4.1[ ) may be evaluated by deforming the contour in the 
complex plane [49J [50] . The deformation is shown in Fig. [6] The left plot (a) shows the Schwarzschild case, and the 
right plot (b) shows the Nariai case. 

On the Schwarzschild spacetime, it is well-known that a 'power-law tail' arises from the frequency integral along 
a branch cut along the (negative) imaginary axis (Fig. |6j part (3)). In the Schwarzschild case, the branch cut is 
necessary due to a branch point in giu,(r, r') at u = [49, 66^ In contrast, for the Nariai case with £ > 0, the 



Wronskian (4.16) is well-defined and non-zero in the limit uj — > 0. For minimal coupling (£ = 0), we find that u = 
is a simple pole of the Green function. In either case, uj = is not a branch point and hence power law decay does 
not arise on the Nariai spacetime. 

The simple poles of the Green function (shown as dots in Fig. [6| occur in the lower half-plane of the complex 



frequency plane. The poles correspond to the zeros of the Wronskian (4.16). The Wronskian is zero when the 
"in" and "up" solutions are linearly-dependent. This occurs at a discrete set of (complex) Quasinormal Mode (QNM) 
frequencies uj q . Beyer [67] has shown that, for the Poschl- Teller potential, the corresponding QNM radial solutions form 
a complete basis at sufficiently late times (t > t c , to be defined below). Completeness means that any wavefunction 
obeying the correct boundary conditions at — > ±oo can be represented as a sum over quasinormal modes, to 
arbitrary precision. Intuitively, we may expect this to mean that, at sufficiently late times, the Green function itself 
can be written as a sum over the residues of the poles. 
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(a) Schwarzschild 



o o o o 



(2) 



o o o o 



(3) 



Re(w) 



(b) Nariai 


Im( w ) 
















\ Re(u>) 












O l 












o ' 










(2) 




\ o o 


o 


o 


o 


o o o / 










— - - 







FIG. 6: Contour Integrals. These plots show the deformation of the integral over frequency to include the poles of the Green 
function (quasinormal modes), for two spacetimes: (a) Schwarzschild [left] and (b) Nariai [right]. 



1. Quasinormal Modes 



Quasinormal Modes (QNMs) are solutions to the radial wave equation (3.11) which are left-going (e iuj p*) at 
p* — > — oo and right-going (e +2u;p *) as p* — > +oo . QNMs occur at discrete complex frequencies uo = uo q for which 

= 0. At QNM frequencies, the "in" and "up" solutions (u 1 ^ and ) are linearly-dependent and the Wronskian 



(4.16) is zero. 

The QNMs of the Schwarzschild black hole have been studied in much detail [68 , 69][70l[7l]. QNM frequencies are 
complex, with the real part corresponding to oscillation frequency, and the (negative) imaginary part corresponding 
to damping rate. QNM frequencies uoi n are labelled by two integers: Z, the angular momentum, and n = 0, 1, . . . oo, 
the overtone number. For every multipole Z, there are an infinite number of overtones. In the asymptotic limit I ^> n, 
the Schwarzschild QNM frequencies approach [52j |72j 175] 



Muo 



(S) 



1 



[±(Z + l/2)-i(n + l/2)] 



(4.19) 



In general, the damping increases with n. The n = ('fundamental') modes are the least-damped. 

The quasinormal modes of the Nariai spacetime are found from the condition = 0. Using (4.13), we find 

1 + (3 — iuoin = —n or — (3 — iuoi n = — n (4.20) 

where n is a non-negative integer. These conditions lead to the QNM frequencies 

ooin = -A-i(n+l/2), (4.21) 

where A is defined in ( |4.5| ). (Note we have chosen the sign of the real part of frequency here for consistency with 
previous studies [521 153] which use a = —u as the frequency variable.) 



2. The Quasinormal Mode Sum 



The quasinormal mode is constructed from (2.6) by t akin g the sum over the residues of th e pol es of gi u (p, p') m 
the complex-cj plane. Applying Leaver's analysis [49 to (2.6), with the radial Green function (4.17), we obtain 



G r Q e ? M (*,/^V;7) = 2Re^ £(2Z + \)Pi{coii)Bi n ui n {p)u ln {pl)e 



■iWln(t — t'—p*—p'^) 



(4.22) 



n=0 1=0 



where the sum is taken over eith er the third or fou rth q uadrant of the frequency plane only. Here, p and t are the 
coordinates in line element (3.8), p* is denned in (3.10), uoi n are the QNM frequencies and B\ n are the excitation 
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factors, defined as 



B ln 



.(out) 



2idi n 



dA 



(in) 



duj 



(4.23) 



and uin(p) are the QNM radial functions, defined by 



UlniP*) 



Xp*) 



(4.24) 



The QNM radial functions are normalised so that ui n (p) —> 1 as p — > 1 (p* — > oo). 

The QNM frequencies have a negative imaginary part; hence the exponentials in (4.22) diverge with n at 'early' 
times t — t' <t c and t c = p* + p*. The exponentials converge with n at 'late' times t — t > t c . Beyer [67] has shown 
that, for late times t — t'> t c , QNMs form a complete basis. Physically, for the QNM sum to be appropriate, sufficient 
coordinate time must elapse for a light ray to propagate inwards from //, reflect off the potential barrier near p = 0, 
and propagate outwards again to p. 



The excitation factors, defined in (4.23), may be shown to be 



Bm = 



1 T(n 



2iui n ) 



2n! [r(l-^ Zn )] 2 ' 

1 r(-n + 2zA) 
2ra! [r(-ra + l/2 + zA)p 



(4.25) 



The steps in the deriva tion a re given in Appendix A. 3 of [53] . 

The Green function (4.22) now takes the form of a double infinite sum, taken over both angular momentum I and 
overtone number n. The convergence of the sum over I is by no means guaranteed. For example, the magnitude of 
Bi n is proportional to (I + l/2) n-1 / 2 in the large-/ limit, for fixed n. Hence, the magnitude of each term in the series 
grows with I. Nevertheless, we will show that well-defined and meaningful values can be extracted from the sums. 



3. Green Function at Arbitrary 'Radii' 



Unfortunately, it is not straightforward to perform the sum o ver n expl icitly for general values of p and p' . Instead 
we must inclu de t he QNM radial functions ui n {p*), defined by (4.7) and (4.24). We n ote that the Legendre function 
appearing in (4.7) can be expressed in terms of a hyper geometric function (see Eq. (4.6)), and the hypergeometric 
functio n can be written as power series about p = 1 (Eq. |4.11 ). At quasinormal mode frequencies, the second term 
i n E q. ( 4.11 ) is z ero, a nd hence the wavefunction is purely outgoing at infinity, as expected. Combining results (4.6), 
(4.7), (4.11) and (4.24) we find the normalised wavefunctions to be 



ui n (p) 



(4.26) 



where p is the radial coordinate in line element (|3.8[) and Si n is a finite series with n terms, 



1 , 1 - p 

-n + 2zA, — n; — n + - + i\; — - — 



*n{p) = t M~". + ?y:" )fc ( = 2F, 



fc! (-n + l/2 + iA) fc 



(4.27) 



where we have adopted the sign convention uoi Tl 
'radii' may be written as the double sum, 



A — i(n + 1/2) of (4.21). Hence the Green function at arbitrary 



G^ M (x,x') = 2Re^(2? + l)P i (cos 7 )e a [ r - 1 «( 2 /(i+P))-in(2/(i + p'))] 

1 71+1/2 



i=0 



n=0 



4e" 



Sln{p)Sln{p') 



(4.28) 
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4- Green Function near Spatial Infinity, p*,p* — > +00 



In the limit that both radial coordinates p* and p^ tend to infinity (p, p r —> 1), the QNM sum (4.22) may be 
rewritten 



lim G^ M (T, 7 ) =2ReY / (2l + l)P l (cos 1 )e iXT f]B ln e-^ +1 



/2)T 



Z=0 



n=0 



where 



T = t — t' — p* — p'* 



(4.29) 



(4.30) 



and A was defined in (4.5). Using the expression for the excitation factors (4.25), the sum over n can be evaluated 
explicitly, as follows, 



1/2 



5>»*- (n+1/2)T = VE 



n=0 



n=0 



T(-n + 2iA) z n 
[r(-n + l/2 + iA)] 2 ^ 



,1/2 



r(2iA) 



2 [T(l/2 + *A)]2§[(l/2 + tA)_ n ]2n! 



(2iA)_ 



(4.31) 



n=0 



where z = e T and (•)& is the Pochhammer symbol. Using the identity (x)- n = ( — l) n /(l — x) n and the duplication 
formula T(z)T(z + 1/2) = 2 1 - 22 v ^r(2z) we find 

00 p -T/2 o2iAp/ 7 -\\ 

^^„e-(" +1 /2)T = f_J_JM_ 2Fl( _ / 3 > _ / 3._ 2) g._ e -T ) (432) 



n=0 



4^ r(l/2 + zA) 



where 2^1 is the hypergeometric function, and j3 was defined in (4.4). Hence the Green function near spatial infinity 
(p, // 1) is 

lim G2T(r >7 ) = ^Reg (* + ^ (4.33) 



5. Green Function Approximation from Fundamental Modes 



Expression (4.28) is complicated and difficult to analyse as it involves a double sum. It would be useful to have 



a simple approximate expression, with only a single sum, which captures the essence of the physics. At late times, 
we might expect that the Green function is dominated by the least-damped modes, that is, the n = fundamental 
quasinormal modes. If we discard the higher modes n > 0, we are left with an approximation to the Green function 
which does indeed seem to capture the essential features and singularity structure. However, as we show in section 



VP it does not correctly predict the singularity times. 
The n = approximation to the Green function is 



^ret X ) 



7r(l + p)(l + p') 



1/2 00 



E ^ + wwi C ° S ^ r (iA) e^+'"(i WV)]^ )(W) , 



1=0 



r(l/2 + iA) 



(4.34) 



Note that the sum over I is approximately periodic in T (exactly periodic in the case £ = 1/8), with period 47r. 



C. Quasilocal Green Function: Hadamard-WKB Expansion 



We now consider the Green function in the quasilocal region, which is needed for the calculation of the quasilocal 
contribution to the scalar self-force given in Eq. (2.4). When spacetime points x and x' are sufficiently close 

together (within a convex normal neighbourhood), the retarded Green function may be expressed in the Hadamard 
parametrix [4T| 142]. 



G ret Or, x) = 6L (x, x) {U (x, x) S (a (x, x')) - V (x, x) (-a (x, x'))} , 



(4.35) 
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where #_ (x, x') is analogous to the Heaviside step- function (unity when x' is in the causal past of x, zero otherwise), 
5(a(x,x')) is the standard Dirac delta function, U (x,x f ) and V (x,x') are symmetric bi-scalars having the benefit 
that they are regular for x' — > x, and a (x, x') is the Synge p~6j [74] [75] world function. Clearly, the term involving 



U(x,x'), the 'direct' part, will not contribute to the quasilocal integral in Eq. (2.4) since it has support only on the 
light-cone, while the integral is internal to the light-cone. We will therefore only concern ourselves with the calculation 
of the function V(x,x'), the 'tail' part, which has support inside the light-cone. 

The fact that x and x' are close together suggests that an expansion of V(x,x') in powers of the separation of the 
points may give a good representation of the function within the quasilocal region. In Ref. [76] we use a WKB method 
(based on that of Refs. [43] [77] [78]) to derive such a coordinate expansion and we also give estimates of its range of 
validity. Referring to the results therein, we have V(x,x') as a power series in (t — t') and (cos 7 — 1), 



+00 

N 2i 



V (x, x') = "a (P) (* - 0* ( cos 7 - 1) j , (4.36) 

i,j=0 

where 7 is the angular separation of the points. In general this expression also includes a third index, fc, corresponding 
to the fc-th power of the radial separation of the points, (p — p') k . However, for the non-radial motion considered in 
the present work, we will only need the terms of order O [(p — p') ] and O [(p — p') 1 ] . The k = terms, Vijo = Vij, 



are given by Eq. (4.36) and Ref. [76] and the k = 1 terms, i^i, are easily calculated from the k = terms using the 
identity [47] 

Viji(p) = -^VijoAP)- ( 4 - 37 ) 



Eq. (4.36) therefore gives the quasilocal contribution to the retarded Green function as required in the present context. 



V. SINGULAR STRUCTURE OF THE GREEN FUNCTION 



In this section we investigate the singular structure of the Green function. We note that one expects the Green 
function to be singular when its two argument points are connected by a null geodesic, on account of the 'Propagation 
of Singularities' theorems of Duistermaat and Hormander [79 , 80 and their application to the Hadamard elementary 
fun ction (which is, except for a constant factor, the imaginary part of the Feynman propagator defined below in 
Eq.( 5.36| )) for the Klein-Gordon equation by, e.g., Kay, Radzikowski and Wald [81 : "if such a distributional bisolution 
is singular for sufficiently nearby pairs of points on a given null geodesic, then it will necessarily remain singular for 
all pairs of points o n tha t null geodesic." 

W e begi n in Sec. V A by exploring the large-/ asymptotics of the quasinormal mode sum expressions (4.33), (4.28) 
and (4.34). The large-/ asymptotics of the mode sums are responsible for the singularities in the Green function. 
We argue that the Green function is singular whenever a 'coherent phase' condition is satisfied. The coherent phase 
condition is applied to find the times at which the Green function is singular. We show that the 'singularity times' 
are exactly those predicted by the geodesic analysis of Sec. Ill C In Sec. VB we introduce two methods for turning 
the sum over / into an integral. We show in Sec. V C that the Watson transform can be applied to extract meaningful 
values from the QNM sums, away from singularities. We show in Sec. VP that the Poisson sum formula may be 
applied to study the behaviour of the Green function near the singularities. We show that there is a four- fold repeating 
pattern in the singular structure of the Green function, and use uniform asymptotics to improve our estimates. In 
Sec. |V E| we rederive the same effects by computing the Van Vleck determinant along orbiting geodesies, to find the 
'direct' part of the Green function arising from the Hadamard form. The two approaches are shown to be consistent. 



A. Singularities of the Green function: Large-1 Asymptotics 

We expect the Green function G ret (x,x') to be 'singular' if the spacetime points x and x' are connected by a null 
geodesic. By 'singular' we mean that G re t(x,x') does not take a finite value, although it may be well-defined in a 
distributional sense. Here we show that the Green function is 'singular' in this sense if the large-/ asymptotics of the 
terms in the sum over / satisfy a coherent phase condition. 
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1. Near spatial infinity p, p' — > +00 



Insight into the occurrence of singularities in the Green function may be obtained by examining the large-/ asymp- 
totics of the terms in the series (4.33). Let us write 



lim G^ m {x,x')=ReY,Qi{T,l) 

1=0 



(5.1) 



The asymptotic behaviour of the gamma function ratio is straightforward: r(iA)/r(l/2 + iX) ~ A x l 2 e 27r / 4 , A 
+00. The large-/ asymptotics of the hypergeometric function are explored in Appendix [Aj We find 



2 F 1 (-f3,-f3;-2f3;-e- T )~ {l + e~ T ) 1/4 exp *A 



In- 



1 



21n2-T 



A — > +00. 



e- 1 -1 

For simplicity, let us consider the special case of spatial coincidence 7 = (near spatial infinity p, p r — > 1), 

-T/2 



Gi(T, 1 = 0) 



(I + 1/2) 



yftil + e-T) 1 /* A 1 / 2 



exp iAln 



1 



1 



27T/4 , A — > 00. 



(5.2) 



(5.3) 



Asymptotically, the magnitude of the terms in this series grows with (/ + 1/2) x / 2 . Hence the series is not absolutely 
convergent. Nevertheless, due to the oscillatory nature of the series, well-defined values can be extracted (see Sec. V C ), 
provided that the coherent phase condition, 



lim arg(ft+i/ft) = 27rJV, N e Z, 



(5.4) 



is not satisfied. In other words, the Green function is 'singular' in our sense if Eq. (5.4) is satisfied. In this case, 

Vl + e- T + l X 



In 



1 



2tt7V, 



N e Z. 



Rearranging, we see that the Green function (4.33) with 7 = is 'singular' if 

T = t - t' - - pi = In [sinh 2 (7rA0] • 



(5.5) 



(5.6) 



Note that the coherent phase condition (5.4) implies that the Green function is singular at precisely the null geodesic 



times (3.23) (with H = 1 and = 27r7V77derived in Sec. VD 



For the more general case where the spacetime points x, x are separated by an angle 7 on the sphere, it is 
straightforward to use the asymptotics of the Legendre polynomials to show that the Green function (4.33) is singular 
when 



T = t - t' - - ^ = In [sinh 2 (A(/>/2)] , where A</> = 2tt7V ± 7 



(5.7) 



again in concordance with (3.23). 



2. 'Fundamental mode ' n — 



In Section IV B 5 we suggested that a reasonable approxima tion t o the Green function may be found by neglecting 
the higher overtones n > 0. The 'fundamental mode' series (4.34) also has singularities arising from the coherent 



phase condition (5.4), but they occur at slightly different times; we find that these times are 



t(n=0) 



A0-ln[(l + p)(l + //)] where = 2nN ± 7. 



(5.8) 



Note that the singularity times T( n=0 ) are periodic. Towards spati al infi nity, (5.8) simplifies to T( n = ) = — 2 In 2, 
which should be compared with the 'null geodesic time' given in (3.23). Clearly, the periodic tim es TL=o) are not 
quite equal to the geodesic times. Neve rthele ss, the latter approaches the former as N — > 00. I n Sec. |VD| we compare 



the singularities of the approximation (4.34) with the singularities of the exact solution (4.33) at spatial infinity 



To investigate the form of the Green function near, and away from, the null cone, we now introduce two methods 
for converting a sum over I into an integral in the complex /-plane. 



18 



Im(i/) 



Re(u) 



FIG. 7: The Watson Transform. The plot shows the contour C\ that defines the Watson transform in Eq. (5.11). Provided the 
integrand is convergent in both quadrants / and IV , the contour may be deformed onto C2 . 



B. Watson Transform and Poisson Sum 



In Sec. IV B[ the distant-past Green function was expressed via a sum over I of the form 



+00 

i=0 



|)Pi(cos 7 ). 



(5.9) 



Here, T{1 + \) may be immediately read off from (4.33), (4.28) and (4.34). The so-called Watson Transform [82] 
and Poisson Sum Formula [83 provide two closely-related ways of transforming a sum over I into an integral in the 
complex /-plane. The two methods provide complementary advantages in understanding the sum over /. 

A key element of the Watson transform is that, when extending the Legendre polynomial to non-integer /, the 
function with the appropriate behaviour is Pi( — cos 7). This is obsc ured by the fact that Pi(— cos 7) = (— l) 1 Pi (cos 
when I is an integer. Our first step is then to rewrite the sum (5.9) as 



+00 

£< 

1=0 



^(Z+|)P / (-cos 7 ), 



(5.10) 



where we have also introduced an integer N for later convenience. Using the Watson transform, we may now express 
the sum (5.9) as a contour integral 



I 



(-1) 



N 



2i 



A2N 



Ci 



™^(z/)P,_ 1/2 (-cos 7 ) 



COs(7TZ/) ' 



(5.11) 



where v — I + 1/2. The contour C\ starts just below the real axis at 00 encloses the points ^=^,1 + ^,2 + . . . 
which are poles of the integrand and returns to just above the real axis at 00. The contour C\ is shown in Fig. [7J If 
the integrand is exponentially convergent in both quadrants / and IV, the contour may be deformed in the complex-/ 
plane onto a contour C2 parallel to the imaginary axis (see Fig. [7|. Note that 'Regge' poles are not present inside 
quadrants / and IV in this case. 

To study the asymptotic behaviour of the Green function near singularities it is convenient to use the alternative 
respresentation of the sum obtained by writing 



1 



COs(7T^) 



1=0 

CO 

-2i^e" i7r(2l+1)(i/_1/2) 



- 2 * E 



i7r(2Z+l)(i/-l/2) 



Im(» > 0, 



Im(V) < 0. 



(5.12) 
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Inserting representation (5.12) into (5.11) leads to the Poisson sum formula, 



+oo poo 

1= Y] (-l) s Re / dve^ iav T(y)P v _ xj2 (co^). 



(5.13) 



The Poisson sum formula is applied in Sec. VP to study the form of the singularities. 



C. Watson Transform: Computing the Series 

Let us now show how the Watson transform may be applied to extract well-defined values from sums over /, even 
though the terms in the series are not absolutely convergent as I — > oo. We will illustrate the approach by considering 
the nth-overtone QNM contribution to the Green function Eq. (4.28) for which case 



T{x,x'\v) = Re 



1 



4e" 



(1 + pXl + p/) 



2+1/2 



g 7 . r iA[T-ln(2/(l+p))-ln(2/(l+p / ))] r (~^ + giA) 

[r(-n + l/2 + iA)] 2 



(5.14) 



x 2 F 1 (-n - 2i\ -n, -n + \ + i\\ (1 - p)/2) x 2 F 1 (-n - 2i\ -n, -n + § + zA; (1 - p')/ 2 ) , 

where A was defined in ( [4^5 1 ). We will choose the integer iV so that this contour may be deformed into the complex 
plane (as shown in Fig. Fn) to a contour on which the integral converges more rapidly. First, we note that the 
Legendre function may be be written as the sum of waves propagating clockwise and counterclockwise, P I/ _ 1 / 2 (cos7) = 

2^- ) i/ 2 ( cos 7) + Q^ ) i/ 2 ( cos 7), where 



2? 

P^z) ± -Quiz) 



(5.15) 



and here Q^{z) is a Legendre function of the second kind. The functions Q^ 1 / 2 nave exponential asymptotics in the 
limit z/7 ^> 1, 



SS /2 (cos 7 ) 



1 



2ttu sin 7 



1/2 



(5.16) 



With these asymptotics, and with the large-/ asymptotics (5.2) of the hypergeometric functions, one finds that the 
contour may be rotated to run, for example, along a line Re(F) = c with c a constant between and \ provided that 
we choose 



N = 



' \(T + log ((p + l)(p' + 1)) + 7 )/(2tt)] for QM 

[(T + log ((p + 1)^ + 1)) + 2ir - 7 )/(2tt)] for Q(") 



(5.17) 



where here [x] denotes the greatest integer less than or equal to x. 

We performed the integrals along Re(z^) = \ an d found rapid convergence of the integrals except near the critical 
times defining the jumps in TV given by Eq. (5.17), when the integrands fall to increasingly slowly. An alternative 



method for extracting meaningful values from series which are not absolutely convergent is described in Sec. VI D 



D. The Poisson sum formula: Singularities and Asymptotics 



In this section we study the singularity structure of the Green function by applying the Poisson sum formula (5.13). 
The first step is to group the terms together so that 



N=0 



where 



Re 



+00 



dvT{y)R^{y^ 7) 



(5.18) 



and 



(_l)iv/2 [ Q (-) i/2 (cos 7 )e^ + QlV 1/2 (cos~f)e 



!(+) 



R 



N 



-iNn 



("I) 



(JV+l)/2 



2L ) i/ 2 ( C0S 7)e 



-i(N+l)irv 



N even, 
N odd. 



(5.19) 
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and Q ] ^_ 1 / 2 were defined in Eq. (5.15). We can now use the exponential approximations for Q ] ^_ 1 / 2 given in (5.16) 
to establish 



R 



(-1) 



N/2 



■N 



(27T^sin7) 1 /2 W_ 1 ^(iV+l)/2 e *7r/4 e ii/((iV+l)7r-7) + cc 



A 7 " even, 
N odd. 



(5.20) 



It should be borne in mind that the exponential approximations (5.16) are valid in the limit 7^ ^> 1. Hence the 
approximations are not suitable in the limit 7^0 case. Below, we use alternative asymptotics (5.28) to investigate 
this case. 



1. 'Fundamental mode' n — 



Let us apply the method to the 'fundamental mode' QNM series (4.34). In this case we have 

1/2 



4e" 



tt(1 + p)(l + p') J r(l/2 + zA) 



o iX X 



where T was defined in (4.30), A was defined in ( |4.5[ ), and 

X = r + ln[(l + p)(l + p') 



Taking the asymptotic limit v — > 00 we find 



4z/e" 



Z7r(l + p)(l + p0 



1/2 



V — > OO. 



Now let us combine this with the 'exponential approximations' (|5.20|) for i?/v 

2e- T 



T{y)R N {y,i) 



1/2 f ( _ 1} iv/2 
7T 2 sin 7(1 + pKl + p') ) I (_i)(^v+i)/2 



e ^(x+(iV+l)7r- 7 ) _ ^ e i^(x-(iV+l)7T+ 7 ) 



for ^ ^ 00. It is clear that the integral in (5.18) will be singular if the phase factor in either term in (5.24) is zero. 



(5.21) 



(5.22) 



(5.23) 



N even, 

TV odd, 

(5.24) 



In other words, each wave Rn gives rise to two singularities, occurring at particula r 'sing ularity times'. We are only 
interested in the singularities for T > 0; hence we may neglect the former terms in (5.24). Now let us note that 



r 

lim / 

e-0+ Jo 



0+ \( + ie 



= i/C + ™*(0 



(5.25) 



Upon substituting (5.24) into (5.18) and performing the integral, we find 

1/2 



r(0) 
N 



2e~ 



sin7(l + p)(l + pf) 



(_l)(^+D/2 



where iff is In in (|5.18|) with T(y) given by (|5.2l|), and 



AO) 



p # +^-ln((l+p)(l+p')) + 



TVtt + 7, 
(JV+1)tt-7, 



N even, 
N odd, 



N even, 
TV odd. 



(5.26) 



(5.27) 



These times t$ are equivalent to the 'per iodic' times identified in Sec. |V A| (Eq. |5.8[ ) and Sec. 



IIIC 



(Eq. 



3.24). 



Let us consider the implications of Eq. (5.26) carefully. Let us fix the spatial coordinates p, p' and 7 and consider 



variations vat — t' only. Each term X^ corresponds to a particular singularity in the mode sum expression (4.34) for 



the (n = 0)-Green function. The Nth singularity occurs at t — t' = tffi. For times close to t^\ we expect the term X^ 



to give the dominant contribution to the (n = 0)-Green function. Eq. (5.26) suggests that the (n = 0)-Green function 
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has a repeating four-fold singularity structure. The 'shape' of the singularity alternates between a a delta-distribution 
(=L5(£ — t f — t$), N even) and a singularity with antisymmetric 'wings' (=Ll/(t — t f — t$), N odd). 

The Nth wave may be associated with the Nth orbiting null geodesic shown in Fig. [2] Note that 'even N' and 
'odd N' geodesies pass in opposite senses around p = (see, for example, Fig. Now, N has a clear geometrical 
interpretation: it is the number of caustics through which the corresponding geodesic has passed. Caustics are 
points where neighbouring geodesies are focused, and in a spheric ally-s ymmetric spacetime caustics occur whenever 
a geodesic passes through angles Acj) = 7r, 27r, 37r, etc. Equation (5.26) implies that the singularity structure of the 



Green function changes each time the wavefront passes through a caustic [84] . 

More accurate approximations to the singularity structure may be foun d by u sing the uniform asymptotics estab- 
lished by Olver [85] (as an improvement on the 'exponential asymptotics' ( 5.16| )), 



Q 



(±) 



-1/2 



(cos 7) 



sin 7 



1/2 



(5.28) 



where H^^) = Jo(-) ±iYo(-) are Hankel functions of the first (+) and second (— ) kinds. This ap proxim ation (5.28) 
is valid in the large-^ limit for angles in the range < 7 < tt. With these asymptotics, we replace (5.20) with 



R 



■N 



sin 7 



1/2 



(-1)*"* p| +1 Ne 
( _ 1)( iv + i)/2 \ H (-)^y(N + i)^ +ir(+)( I/7 ) e -<(^+i> 



H ( -\^)e 



-iN-KV 



N even, 
TV odd. 



(5.29) 



In Appendix [B] we derive the following asymptotics for the 'fundamental mode' (n = 0) Green function, 

2^(7) 

T + 7 )_ x ]l/2-V- 

2tt 



r(0) 



^[(2.- 7 )-x][(2. + 7T^F i?(27/[(27r + 7) 



z 



(o) 



2[x-(2^- 7 )] 3 / 2 

r 27T 



2 F 1 3/2, 1/2; 2; 



•7 



X - 2tt + 7 



X < 2tt- 
X > 2tt- 



7, 
■7» 



^(7)*(x-(2tt + 7 )), 
^ (7)V ^ 1/2; 2;^ 



27T — 7 



X<2tt- 
X > 2tt- 



l 2[X " (2tt - 7 )] 3 / 2 *" ^ V" X ~ 2tt + 7 

where £^ is the complete elliptic integral of the second kind, x was defined in (5.22) and 

1/2 



"7, 
"7, 



(5.30) 



(5.31) 



7 



sin 7 



Tra+pxi+po 



1/2 



(5.32) 



The asymptotics ( |5.30| ) and (|5.31| ) provid e insight into the singularity structure near the caustic at Acj) = 2tt. Figure 
[8] shows the asymptotics ( |5.30| ) and (5.31) for two cases: (i) 7 = 7r/20 (left) and (ii) 7 = (right). In the left plot, 
the x[ ^ integral has a (nearly) antisymmetric form. The integral is a delta function with a 'tail'. However, the 

'tail' is exactly cancelled by the integral in the regime x > 27r + 7. The cancellation creates a step discontinuity 
in the Green function at x = 27T + 7. The form of the divergence shown in the right plot (7 = 0) may be understood 
by substituting 7 = into (5.30) to obtain 



2<°>(7 = 0) 



(i + p)(i + p') 



1/2 



(2tt-x)- 3/2 . 



(5.33) 



2. Near spatial infinity, p, p' 



It is straightforward to repeat the steps in the abov e analysis for the closed- form Green function (4.33), valid for 
p, p' — > 1. We reach a result of the same form as ( |5.26[ ), but with modified singularity times, 



.(») = , / , fin (smh 2 ([Nir + 7 ]/2)) , N even, 

P * +P * + 1ln(sinh 2 ([(iV+l)7r-7]/2)), N odd, 



(5.34) 
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Singularities of Green function from Poisson Sum Formula Singularities of Green function from Poisson Sum Formula 
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FIG. 8: Singularities of the 'fundamental mode' Green function (4-34]) nea r the caust ic at 2tt. These plots show the X\ (dashed) 
and X2 (dotted) contributions to the Poisson sum, given explicitly by ( |5.3Q ) and ( 5.31 ). The left plot shows an angular separation 
7 = 7r/20 and the right plot shows angular coincidence 7 = 0. Note that, for T > 2tt + j — 2 In 2, the ii and I2 integrals are 
equal and opposite and will exactly cancel out (see text). 



corresponding to the geodesic times (3.23). For instance, with the exponential asymptotics (5.20) applied to (4.33) 
we obtain 



r (oo) 
N 



2 sin 7 + e~ T 



1/2 



\-l) N ^ 2 s(t-t f -t^), N even, 

(_l)(A / "+l)/2 



,(oo) 
N 



odd. 



(5.35) 



In Sec. VII the asymptotic expressions derived here are compared against numerical results from the mode sums. 



We believe that the 4-fold cycle in the singularity structure of the Green function which we have just unearthed 
using tricks we picked up from seismology [83] is characteristic of the S 2 topology (different types of cycle arising in 
different cases). This cycle may thus also appear in the more astrophysically interesting case of the Schwarzschild 
spacetime. Since this cycle does not seem to be widely known in the field of General Relativity (with the notable 
exception of [84]), in Appendix |C| we apply the large-/ asymptotic analysis of this section to the simplest case with 
S 2 topology: the spacetime of T x S 2 , where the same cycle blossoms in a clear manner. 



E. Hadamard Approximation and the Van Vleck Determinant 



In this section, we rederive the singularity structure found in (5.26) and (5.35) using a 'geometrical' argument based 
on the Hadamard form of the Green function. In Sec. |IV C] we used the Hadamard parametrix of the Green function 
to find the quasilocal contribution to the self-force. Strictly speaking, the Hadamard parametrix of Eq. ( 4.35| ), is only 
valid when x and x' are within a convex normal neighborhood [41]. Nevertheless, it is plausible (particularly in light of 
the previous sections) that the Green function near the singularities may be adequately described by a Hadamard-like 
form, but with contributions from all appropriate orbiting geodesies (rather than just the unique timelike geodesic 
joining x and x'). 

We first intr oduc e the Feynman propagator Gf(x^x') (see, e.g., [11, 86J) which satisfies the inhomogeneous scalar 
wave equation (|2.3|). The Hadamard form, which in principle is only valid for points x' within the normal neighbour- 
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hood of x, for the Feynman propagator in 4-D is [TTJ [87] 

U (x, x' 



Gf(x, x') 



% 

2^ 



+ V(x, a:') ln(cr + ie) + a:') 



(5.36) 



IV C) and VF(a;,a/) are bitensors which are regular at coincidence 



where U(x,x f ), V(x,x f ) (already introduced in Sec 
(x — > a; 7 ) and a(x,x') is Synge's world function: half the square of th e geodesic distance along a specific geodesic 
joining x and x' . Note the Feynman prescription L o~ — > a + ie' in ( |5.36| ), where e is an infinitesimally small positive 
value. 

The retarded Green function is readily obtained from the Feynman propagator by 



G ret (x,x') = 26-(x,x') Re (G F (x, x')) , 



(5.37) 



which yields (4.35) inside the convex normal neighbourhood, since U(x,x'), V(x,x f ) and W(x,x') are real-valued 
there. We posit that the 'direct' part of the Green function remains in Hadamard form, 

U(x,x') (&{&) + 





~.u(x,x f y 
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= Re 


lim —Re 


e^0+ 7T 


a + ie 





— I 



7TO~ 



(5.38) 



even outside the convex normal neighbourhood [note that outside the convex normal neighbourhood we still use the 
term 'direct' part to refer to the contribution from the U(x,x') term, even if its support may not be restricted on 
the null cone anymore]. It is plausible th at th e Green function near the iVth singularity (see previous section) is 
dominated by the 'direct' Green function (5.38) calculated along geodesies near the 7Vth orbiting null geode sic. To 
test this assertion, we will calculate the structure and magnitude of the singularities and compare with (5.35). 
In four dimensional spacetimes, the symmetric bitensor U(x,x') is given by 



U(x,x') = A 1/2 (x,a/), 



(5.39) 



where A(x,x / ) is the Van Vleck determinant [88, 89, 90 j. The Van Vleck determinant can be found by integrating a 
system of transport equations along the appropriate geodesic joining x and x' . The first of these [16], 



A- 



dlnA 

~d\~ 



(5.40) 



is a transport equation for the Van Vleck determinant itself, with the initial condition A(A = 0) = 1. Here, A is an 
affine parameter along the geodesic joining x and x' \ and cr a p = V pV a o- is the second covariant derivative (taken with 
respect to spacetime point x) of Synge's world function, which in turn is found from the coupled system of transport 
equations [57} l9T] 



A 



d(j c 



dX 



(5.41) 



and the boundary condition cr a ^(A = 0) = 5%. 

In principle, transport equations ( |5.40 ) and (5.41) may be integrated numerically to determine the Van Vleck 
determinant along any given geodesic on any given spacetime. This is the approach that we might take on, for 
example, in Schwarzschild. A numerical approach is not necessary for the Nariai spacetime, however. This spacetime 
is the Cartesian product of a two-sphere with a 2-D de Sitter spacetime. On a product spacetime Ai = M.\ x M. 2 we 
may make the following decomposition: 



o-i -hcr 2 , 



A = AiA 2 , 



(5.42) 



where <r^ and A^ (i = 1,2) are, respectively, Synge's world function and the Van Vleck determinant on the manifold 
Aii. It will be shown in a forthcoming work [92] that the Van Vleck determinant on the Nariai spacetime when the 
two points x and x' are within the normal neighbourhood is simply 



A(x,x / ) 



7 
sin 7 



V 

sinh?? 



(5.43) 



where 7 G [0, 7r) is the geodesic distance traversed on the two-sphere and 77 is the geodesic distance traversed in the 
two-dimensional de Sitter subspace. Hence the Van Vleck determinant is singular at the angle 7 = tt. 
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This may be seen another way. Using the spherical symmetry, let us assume without loss o f generality that the 
motion is in the plane = const, from which it follows that the equation for a^^ in (5.41) decouples from the 
remainder; 



it is 



dO 



= o 



(5.44) 



Here we have rescaled the affine parameter A to be equal to the angle subtended by the geodesic. Note that here 
we let take values greater than 7r. It is straightforward to show that the solution of Eq. (5.44) is a^^ = OcotO. 
Hence <7^ is singular at the angles 6 = n, 2tt, 37r, etc. In other words, the Van Vleck determinant is singular at 
the antipodal points, where neighbouring geodesies are focused: the caustics. The Van Vleck determinant may be 
separated in the following manner: A = A^A^, where 



dlnA 



ty 



Eq. (5.45) yields 



dO 



In Ad, = In 



f 



d0' cot 6' 



(5.45) 
(5.46) 

(5.47) 



which can be integrated analytically by following a Landau contour in the complex 6^-plane around the (simple) poles 
of the integrand (located at 0' = kir, k G Z), which are the caustic points. Following the Feynman prescription 
'<j — > a + ie\ we choose the Landau contour so that the poles lie below the contour. We then obtain (setting # = 
without loss of generality) 



A,* = 



sm( 



(5.48) 



Here, N is the number of caustic points the geodesic has passed through. The phase factor, obtained by continuing 
the co ntour of integration past the singularities at = 7r, 27r, etc., is crucial. Inserting the phase factor e~ lN 7T / 2 in 
( |5.38 ) leads to exactly the four-fold singularity structure predicted by the large-/ asymptotics of the mode sum (5.35). 
That is: 
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(5.49) 



7T<7 



The accumulation of a phase of V on passing through a caustic, and the alternating singularity structure which 
results, is well-known to researchers in other fields involving wave propagation - for example, in acoustics [93] , 
seismology [83] , symplectic geometry [94] and quantum mechanics [95] the integer N is known as the Maslov index [96j 
[97]. 

We would expect to find an analogous effect in, for example, the Schwarzschild spacetime. The four-fold structure 
has been noted before by at least one researcher [84]. Nevertheless, the effect of caustics on wave propagation in 
four-dimensional spacetimes does not seem to have received much attention in the gravitational literature (see [98] [99] 

for exceptions). 

To compare the singularities in the mode-sum expression ( |5.35 ) with the singularities in the Hadamard form (5.49), 
let us consider the 'odd-n' singularities of 1/a form. We will rearrange (5.35) into analogous form by expanding a to 
first-order in t — t^°^ , where t^°^ 
For th e orb iting geodesies described in Sec. 
using (3.23) yields 



is the Nth singularity time for orbiting geodesies starting and finishing at p 

H h2 ■ 



mc 



we have a ■ 



1. 

1)0 2 . At p — > 1, expanding to first order and 
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-(HO) tanh(#<9/2) (t - ^° } ) . (5.50) 
The mode-sum expression (5.35) may then be rewritten in analogous form to the L N odd' expression in ( 5.49| ), 
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Hr]s'mh(r]/2)\ 1/2 
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(5.51) 



Here r] = HQ, where H is the constant of motion introduced in Sec. Ill C We find very good agreement between 



(5.43) and (5.51) in the > tt regime. The disagreement at small angles is not unexpected as the QNM sum is invalid 



at early times (or equivalently, for orbiting geodesies which have passed through small angles 6). 
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VI. SELF-FORCE ON THE STATIC PARTICLE 



In this section we turn our attention to a simple case: the self-force acting on a static scalar particle in the Nariai 
spacetime. By 'static' we mean a particle with constant spatial coordinates. It is not necessarily at rest, since its 
worldline may not be a geodesic, and it may require an external force to keep it static. In Sec. | VI A] we review previous 
calculations for the static self- force on a range of spacetimes and in Sec. VI B we explore one such analytic method for 
computing the static self-force in Nariai spacetime. This method, based on the massive field approach of Rosenthal 
[40] provides an independent check on the matched-expansion approach. In Sec. VIC we describe how the method 
of matched expansions may be applied to the static case. To co mpute the self-force, we require robust nume rical 
methods for evaluating the quasinormal mode sums such as (4.28); two such methods are outlined in Sec. VI D The 



results of all methods are validated and compared in Sec. VII 



A. The Static Particle 



A static particle - a particle with constant spatial coordinates - has been the focus of several scalar self-force 
calculations, in particular for the Schwarzschild spacetime [37j |40j [47l I1QQ[ \101\ 1102} I103j . Although it may not be 
a particularly physical case, it is frequently chosen because it involves relatively straightforward calculations and has 
an exact solution for the Schwarzschild spacetime. It therefore provides a good testing ground for new approaches to 
the calculation of the self-force. 

Smith and Will lOT] calculated the self-force on a static electric charge in the Schwarzschild background and found 
it to be non-zero. In |100j , Wiseman considered the analogous case of a static scalar charge in the case of minimal- 
coupling (i.e., £ = 0) in Schwarzschild. Using isotropic coordinates, he managed to sum the Hadamard series for the 
Green function in the static case (i.e., the "Helmholtz"-like equation in Schwarzschild with zero-frequency, u = 0) 
and thus obtain in closed form the field created by the static charge in the scalar and also electrostatic (already found 
in |104| 1105] using a different method) cases. He then found the self-force to be zero in the scalar, minimally-coupled 
case. 

In [103] , the calculation of the self- force on a static scalar charge in Schwarzschild is extended to the case of 
non-minimal coupling (£ ^ 0) and is found to be zero as well. The fact that the value of the scalar self-force in 
Schwarzschild is the same (zero) independently of the value of the coupling constant is in agreement with the Quinn- 
Wald axioms [14] [15]: their method relies only on the field equations, and these are independent of the coupling 
constant in a Ricci-flat spacetime such as Schwarzschild. The calculation (without using the Quinn-Wald axioms) is 
by no means trivial, however, since the effect of the coupling constant might be felt through the stress-energy tensor 
(in fact, [103] corrected a previous result in [106, 11107] . where the self- force had been incorrectly found to be non-zero). 
Rosenthal [40] has also considered this case of a static particle in Schwarzschild and used it as an example application 
of the massive field approach [108 to self-force calculations. 

On the other hand, Hobbs |109] showed that, in a conj vrmally- flat spacetime, the "tail" contribution to the self-force 
on an electric charge (on any motion, static or not) is zero. The only possible contribution to the self-force might 
then come from the local Ricci-terms, which are zero in cases of physical interest such as in de Sitter universe. 

Noting the conformal-invariance of Maxwell's equations, one would then expect the "tail" contribution to the scalar 
self-force to also be zero in the two following cases: (1) for a charge undergoing any motion in a conformally-flat 4-D 
spacetime with conformal-coupling (i.e., £ = 1/6), and (2) for a static charge (where the time-independence effectively 
reduces the problem to a 3-D spatial one), in a spacetime such that its 3-D spatial section is conformally-flat and 
with conformal-coupling in 3-D (i.e., £ = 1/8). Indeed, in a recent article [110] it was shown that the scalar self-force 
on a massless static particle in a wormhole spacetime (with non-zero Ricci scalar and where the 3-D spatial section is 
conformally-flat) changes sign at £ = 1/8 and it is equal to zero at this 3-D conformal value. 

The Nariai spacetime, not being Ricci-flat and being conformal to a wormhole spacetime (and so with conformally- 
flat 3-D spatial section), suggests a very interesting playground for calculating the self-force: What role does the 
coupling constant £ play? Do particular values such as £ = 1/6 (4-D conformal-coupling) and £ = 1/8 (3-D conformal- 
coupling, so a particular value in the case of a static charge) yield particular values for the self-force? Do they support 
the Quinn-Wald axioms? 



B. Static Green Function Approach 

The conventional approach to calculating the self-force on a static particle due to Wiseman [100] uses the 
'scalarstatic' Green function. Following Copson [T04] , Wiseman was able to obtain this Green function by sum- 
ming the Hadamard series. Only by performing the full sum was he able to verify that his Green function satisfied 
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the appropriate boundary conditions. Linet [111] has classified all spacetimes in which the scalarstatic equation is 
solvable by the Copson ansatz and the Nariai metric does not fall into any of the classes given. Therefore, instead 
we work with the mode form for the static Green function. This corresponds to the integrand at uj = of Eq. (4.18) 
(with integral measure |^), 



G static^ ft; /A ft') = ]T(2Z + l)fl(cos 7 ) 



1=0 



nP-i/2+i\(-p<)P-i/2+i\(p>) 
2 cosh(7rA) 



(6.1) 



where, as before, A = y(/ + ^) 2 + d. This equation, having only one infinite series, is amenable to numerical 
computation. 

To regularise the self- force we follow the method of Rosenthal [40] , who used a massive field approach to calculate 
the static self-force in Schwarzschild. Following his prescription, we calculate the derivative of the scalar field and of 
a massive scalar field. In the limit of the field mass going to infinity, we obtain the derivative of the radiative field 
which is regular. This method can be carried through to Nariai spacetime, where it yields the expression 



ma p = q 2 (l - p 2 )* lim 



dpG sta tic(p, ft; p' \ ft) 



(p-p') 2 i-p 2 



The singular subtraction term may be expressed in a convenient form using the identity (112} 1113] , 
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calculated sums/ integral 



Subtracting ( |6.3[ ) from ( |6.1[ ), we can express the regularised Green function as a sum of two well-defined and easily 

G s tatic(p, ft; p', ft) 1 = Z(p, p') + J(p, p') 
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(6.5) 



(6.6) 



J(p, p r ) may either be calculated directly as a sum or by using the Watson-Sommerfeld transform to write 



(l + i j = $£e - f dz tan(7rz)g(z) 



where 7 runs from to 00 just above the real axis, and for us 
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The contribution from the second term can be best evaluated by rotating the original contour to a contour 7', running 
from to zoo just to the right of the imaginary axis. This is permitted since the Legendre functions are analytic 
functions of their parameter and the contribution from the arc at infinity vanishes for our choice of g(z). From the 
form of g(z) it is clear that it possesses poles along the contour 7' but these give a purely imaginary contribution to 
the integral. We conclude that 
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(6.11) 



where V denotes the Principal Value. These integrals and that defining T(p,p') and their derivatives with respect to 
p are very rapidly convergent and easily calculated. 



C. Matched Expansions for Static Particle 



In Sec. [TTJ we outlined the method of matched expansions. In this subsection, we show how to apply the method 
to a specific case: the computation of the self-force on a static particle in the Nariai spacetime. 
The four velocity of the static particle is simply 



0, ^ = (l-p 2 ) 



-1/2 



and hence dr' = (1 - p 2 ) x l 2 dt' . We find from Eqs. (1.4), (1.5) and (2.1) that ma 1 = ma 



and 
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where denotes p artial differentiation with respect to the coordinate x^. We note that in the tail integral of the 
mass loss equation, (6.14), the time derivative d t may be replaced with — d t > since the retarded Green function is a 



function of (t — t r ). Hence we obtain a total integral, 
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The total integral depends only on the values of the Green function at the present time and in the infinite past 
(t f — > 00). The QNM sum expressions for the Green function (e.g. Eq. (4.33)) are zero in the infinite past, as the 
quasinormal modes decay exponentially. The value of the Green fu nction at coincidence (t f — » t) is found from the 
coincidence limit of the function — V(x,x') in the Hadamard form (4.35). It is ^(1 — 6£)i?, which exactly cancels 



the local contribution in the mass loss equation (6.14). It is no surprise to find that this cancellation occurs - the 
local terms were originally derived from the coincidence limit of the Green fu nctio n. In fact, because = due 
to time-translation invariance, we can see directly from the original equation (1.5) that the mass loss is zero in the 
static case. 



Now let us consider the radial acceleration (6.13). The acceleration keeping the particle in a static position is 



constant (a p = 0). The remaining tail integral may be split into two parts, 
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For the first part of ( |6.16| ), we use the quasilocal calculation of V(x,x') from Sec. IV C As V(x,x') is given as a 
power series in (p — p') and (t — £'), the derivatives and integrals can be done termwise and are straightforward. The 
quasilocal integral contribution is therefore simply 
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The second part of (6.16) can be computed using the QNM sum (4.28). To illustrate the approach, let us rewrite 
p~28} as 



Uln(p)- 



(6.18) 



Applying the derivative with respect to p and taking the integral with respect to t' leads to 
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(6.19) 



It is straightforward to find the derivative of the radial wavefuncti on fr om the definition (4.26). In Sec. 
outline two methods for numerically computing mode sums such as (6.19). 



VI D 



The self- force comput ed via (6.16), (6.17) and (6.19) should be independent of the choice of the matching time (we 
verify this in Sec. VII D ). This invariance provides a useful test of the validity of our matched expansions. Additionally, 



through varying At we may estimate the numerical error in the self- force result. 



D. Numerical Methods for Computing Mode Sums 



The static-self-force calculation requires the numerical calculation of mode sums like ( |6.19| ). We used two me thod s 
for robust numerical calculations: (1) 'smoothed sum', and (2) Watson transform (described previously in Sec. VC). 
We see in Sec. IVIII that the results of the two methods are consistent. 

The 'smoothed sum' method is straightforward to describe and implement. Let us suppose that we wish to extract 
a numerical value from an infinite sum 



J2 ai 

1=0 

which may not be absolutely convergent (i.e. |a/ + i/a/| > 1). We may instead compute the finite sum 



(6.20) 
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(6.21) 
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where is large enough to suppr ess any high-/ oscillations in the result (typically > 4/ cut ). We find that (6.21) 
is a good approximation to (6.20) provided we are not within St ~ 1/Z cut of a singularity of the Green function. 
Increasing the cutoff / cut therefore improves the resolution of the singularities. 



VII. RESULTS 



We now present a selection of results from our numerical calculations. In Sec. VII A the distant past Green function 
is examined. We plot the Green function as a function of coordinate time t — t' for fixed spatial points. A four-fold 
singularity structure is observed. In Sec. |VIIB we test the asymptotic approximations of the singular structure, 
derived in Sees. VE and VD (Eqs. 5.35 and 5.49). We show that the 'fundamental mode' (n = 0) series (4.34) is 



a good approximation of the exact result (4.33), if a 'time-offset' correction is applied. In Sec. VII C the quasilocal 



and distant past expansions for the Green function are compared and matched. We show that the two methods for 
finding the Green function are in excellent agreement for a range of matching times Ar. In Sec. |VII D| we consider 
the special case of the static particle. We present the Green function, the radiative field and the self-force in turn. 
The r adial s elf-force acting on the static particle is computed via the matched expansion method (de scribe d in Sees. 
[D] and [yiCl ) , and plotted as a function of coordinate p, and compared with the result derived in Sec. VI B 
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Y = o 

£ = 1/6 



8 10 

T = t-t'-p.-p.' 



FIG. 9: Distant Past Green Func tion for spatially- coincident points near infinity (p — p' — > 1, 7 = 0). The Green function was 
calculated from mode sum (4.33) numerically using the smoothed sum method (6.21) with Z cut = 200 and curvature coupling 
factor £ = 1/6. 



A. The Green Function Near Infinity from Quasinormal Mode Sums 



Let us begin by looking at the Green function for fixed points near spatial infinity, p = p r — > 1 (i.e. p* 



+00). The Green function ma y be c omputed numerica lly by applying either the Watson transform (Sec. |VC ) or the 
'smoothed sum' method (Sec. VI D) to the QNM sum (4.33). 

Figure 9] shows the Green function for f ixed s patially-coincident points near infinity (p = p' — > 1, 7 = 0). The Green 
function has been calculated from series (4.33) using the 'smoothed sum' me thod . It is plotted as a function of QNM 
- t' — (p* + pi). We see that singularities occur at the times (3.23) predicted by the geodesic analysis 
In this case, Tq = ln[sinh 2 (A/~7r)] 



t 



time, T 
of Sec. 
T w 4.893, t 



IIIC 



4.893, 11.180, 17.463, etc. At times prior to the first singularity at 
le Green function shows a smooth power-law rise. At the singularity itself, there is a feature resembling 
a delta-distribution, with a negative sign. Immediately after the singularity the Green function falls close to zero 
(although there does appear a small 'tail'). This behaviour is even more marked in the case £ = 1/8 (not shown). A 
similar pattern is found close to the second singularity at T « 11.180, but here the Green function takes the opposite 
sign, and its amplitude is smaller. 

Figure 10 shows the Green function for points near infinity (p = p' — > 1) separ ated b y an angle of 7 = tt/2. The Green 
function shown here is computed from the 'fundamental mode' approximation (4.34), again using the 'smoothed sum' 

etc, 



method. In this case, the singularities occur at 'periodic' times ( |3.24| ), given by T 
where Ac/) 



7r/2, 37r/2, 57r/2, . . .. As discussed, there is a one-to-one corres 



orbiting null geodesies, and the four-fold singularity pattern predicted in Sec. VD 



2 In 2 « 0.1845,3.326,6.468, 
p onde nce betw een singula rities and 
KM and Sec. [VE| ([5735]) is clearly 



visible. Every 'even' singularity takes the form of a delta distribution. Numerically, the delta distribution is manifest 
as a Gaussian-like spike whose width (height) decreases (increases) as Z cut is increased. By contrast (for 7 7^ 0,7r), 
every 'odd' singularity diverges as 1/(T — T c ); it has antisymmetric wings on either side. The singularity amplitude 
diminishes as T increases. 



B. Asymptotics and Singular Structure 



The analyses of Sec. 



V D| and Sec. |VE| yielded approximations for the singularity structure of the Green function. 



In particular, Eq. (|5.35|) gives an estimate for the amplitude of the 'odd' singularities as p = p' — > 00. We tested our 

Figure [TT| shows the Green function near the singularity associated 
37r/2. The left plot compares the numerically-determined Green 



numerical computations against these predictions, 
with the null ge odesic passing through an angle A(f) ■ 
function (4.34) with t he asy mptotic prediction ( 5.35[ ). The right plot shows the same data on a log- log plot. The 
asymptotic prediction (5.35) is a straight line with gradient —1, and it is clear that the numerical data is in excellent 
agreement. 
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FIG. 10: Distant Past Green Function near spatial infinity (p — p' — » 1 ) for points separated by angle 7 = 7r/2. The Green 
function was calculated from the 'fundamental mode' approximation (4.34) numerically using the smoothed sum method with 
lent = 1000. Note the four- fold singularity structure (see text). 



Mode sum l max = 2500 
Approximation 



10 - { - 2 - 




FIG. 11: Green function near the singularity arising fro m a nu ll geodesic passing through an angle Acj) = 3tt/2 an d with p = 
p — >• oo. The 'fundamental mode' (n — 0) Green function (4.34) (with Z C ut = 2500) is compared with approximations (5.35) and 
( 5.49| ) from considering high-/ asymptotics. The approximations give G re t ~ -0.04268/(T - T c ) and G re t ~ -0.04344/ (T - T c ), 
respectively. The left panel shows the Green function in the vicinity of the ('periodic') singularity at T c = 3tt/2 — 2 In 2 w 3.3261. 
The right panel shows the same data on a log- log scale, and compares the mode sum (dashed) with the approximation (dotted). 
The discrepancy close to the singularity may be improved by increasing Z cut . 



Improved asymptotic expressions for the singular structure of the fundamental mode Green function were given in 
and ( 5.31 ). These asymptotics are valid all the way up to 7 = 0. Figure 12 compares the asymp totic expressions 

It is clear that 



5.30 



5.30 



and (5.31 ) (solid line) with numerical computations (broken lines) from the mode sum (4.34 



the asymptotics (5.30) and (5.31) are in excellent agreement with the numerically-determined Green function. Closest 
agreement is found near the singular times, but the asymptot ics pr ovide a remarkably good fit over a range of t. 

) 'fundamental mode' (n = 0) approximation (4.34) is compared with the exact QNM Green function 
Away from singularities, the former is found to be a good approximation to the latter. However, close to 



In F ig. 13 



singularities this is not the case. The singularities of the 'fundamental mode' approximation (4.34) occur at slightly 
differe nt tim es to the singularities of the exact soluti on (14.33 ), as discussed in Sec. V A| For t he fun damental mode 
series (4.34), the singularity times T^ eg given in Eq. (5.8) are periodic. For the exact solution (4.33), the singularity 
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-1.5 



FIG. 12: S ingula rities of the 'Fundamental Mode' Green function \4 . 3J$ compared with asymptotics from the Poisson sum 
\5. 30\ and \5.31\ . The left plot shows a small angular separation 7 = 7r/20, and the right plot shows coincidence 7 = 0, for 
p,p -> 1. 



times T^ xact are precisely the 'null geodesic times' given in Eq. (3.23) 



Remarkably, if we apply a singularity time offset 
to the 'fundamental mode' approximation (T — > T + AT where AT = T^ xact — T^ eg ) we find that the 'fundamental 
mode' Green function is an almost perfect match to the exact Green function. This is clearly shown in the lower plot 
of Fig. [l3] Comparing the series ( 4.34[ ) and (4.33) we see that, in both cases, the magnitude of the terms in the series 
increases as (I + 1/2) 1 / 2 in the large-/ limit. This observation raises the possibility that the n = modes may give the 
essential features of the full solution; if true, this would certainly aid the analysis of the Schwarzschild case, where it 
is probably not feasible to perform a sum over n analytically. 



C. Matched Expansions: Quasilocal and Distant Past 



Let us n ow tu rn our attention to the match between quasilocal and distant past Green functions. The quasilocal 
expansion (4.36) is valid within the convergence radius of the series, t — t' < tQL, while the QNM sum is convergent at 
'late' times, t — t' > p* + p^. Hence, a matched expansion method will only be practical if the quasilocal and distant 
past Green functions overlap in an intermediate regime p* + < t — t' < £ql- It is expected that the convergence 
radius of the quasilocal series, £ql, will lie within the normal neighbourhood, t^N, of spacetime point x. The size 
of the normal neighborhood is limited by the earliest time at which spacetime points x and x' may be connected by 
more than one non-spacelike geodesic. Typically this will happen when a null geodesic has orbited once, taking a 
time tNN > P* + P*, so we can be optimistic that an intermediate regime will exist. To test this idea, we computed 



the quasilocal Green function using (4.36), and the distant past Green function (4.28) for a range of situations. 

Figure 14 shows the retarded Green function as a function of coordinate time t—t' for a static particle at p = p' = 0.5. 
At early times, the quasilocal Green function is well-defined, but the distant past Green function is not. Conversely, 
at late times the quasilocal series is not convergent. At intermediate times 1.099 < At < 3.45, we find an excellent 
match. Figure [14] also shows that the results of the two numerical methods for evaluating QNM sums are equivalent. 
That is, the Green function found from the Wats on tran sform (Sec. V C t red line) coincides with the Green function 
calculated by the method of smoothed sums (Sec. VI D black dots). 



Let us now examine the matching procedure in more detail. Figure 15 shows the match between the distant past 
and quasilocal Green functions, computed from ( |4.28 ) and (4.36), in the case p = p' = 0.5. In Fig. 15 , the left plot 
shows the case for conformal coupling £ = 1/6 and the right plot shows the case for £ = 1/8. Note that Green function 
tends to the constant value ^(1 — 6(,)R in the limit At — > + . In both cases, we find that the fit between 'quasilocal' 
and 'distant past' Green functions is good up to nearly the radius of convergence of the quasilocal series. 

Figure [l6| quantifies the accuracy of the match between quasilocal and distant past Green functions. Here, we have 
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FIG. 13: Singularities of 'Fundamental Mode' approximation. 
(p = p — > 1, p* = pi — > oo ), for an angular separation 7 
approximation from the fundamental modes [blue 



The top plot shows the Green function near spatial infinity 
= 7r/2. It compares the exact solution (4.33) [red] and the 
as a function of time T — t — t' — p* — pi . The singularities occur at two 



distinct times, marked T eX act and T reg respectively. If a time offset is applied to the 'Fundamental Mode' approximation (see 
text), then we find the singularities look remarkably similar (lower plot). 



Matching Quasi-local and Non-local Green Functions (£=1/8) 



-- Quasilocal (60th order) 






• Smoothed Sum (n=0) 






— Watson Transform (n=0) 













FIG. 14: Matching of the Quasilocal and Distant Past Green Functions for £ = 1/8 and for a static particle at p = p — 0.5. 
Here, the Green functions are plotted as functions of coordinate time, t — t' (not 'QNM time' as in most other plots). The 
black dots and red line show the results of the 'smoothed sum' and 'Watson transform' methods applied to compute the QNM 
sum (4.34) (see text). The blue dashed line shows the quasilocal series expansion taken to order (t — t') 60 . The distant past 
Green function cannot be computed for early times t — t' < p* + pi = 2tanh _1 (l/2) = 1.099, whereas the quasilocal series 
diverges at large t — t' > 3.46. In the intermediate regime, we find excellent agreement (see also Figs, 
the quasilocal Green function tends to ^(1 — 6^)R — ^ m ^ ne limit t' — > t. 
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and 16). Note that 
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Comparison of QL and DP Green Functions 



Quasilocal series (60th order) 




QNM sum (l cut = 100, n <= 4) 




S = 1/6 




angley = J 




p = p' = 0.5 /] 









Comparison of QL and DP Green Functions 





Quasilocal series (60th order) 

QNM sum (l cut = 100, n <= 4) 


^=1/8 






angle y = 






p = p' = 0.5 







FIG. 15: Matching of the Quasilocal and Distant Past Green Functions. The left plot shows curvature coupling £ = 1/6 and 
the right plot shows £ = 1/8, for a static particle at p — p — 0.5. Note the timescale on the horizontal axis, T — t — t' — p* — p'*. 
The distant past Green function cannot be computed for T < 0, whereas the quasilocal series clearly diverges at large T. In 
the intermediate regime, we find excellent agreement. 




FIG. 16: Error in Matching the Quasilocal and Distant Past Green Functions. This plot shows the difference between the 
quasilocal and distant past Green functions in the matching regime. The magnitude of the Green function as a function of 
time T = t — t' — p* — p'* is shown as a solid red line. The broken line s show the ' ma tching error': the difference between the 
quasilocal and QNM sum Green functions, for various l cut and n max in (6.21) (see Sec. VI D). Note the logarithmic scale on the 
vertical axis. It is clear that the matching error is reduced by increasing n max and Z C ut, and that the best agreement is found 
close to the radius of convergence of the quasilocal series (at T ~ 1.9). The plot shows that matching accuracy of above one 
part in 10 4 is achievable. 



used the 'smoothed sum' method (Sec. VI D) to compute the distant past Green function from (4.28). To apply this 
method, we must choose appropriate upper limits for / (angular momentum) and n (overtone number). We have 
experimented with various cutoffs / cut and n max . As expected, better accuracy is obtained by increasing / cut and 
n max , although the run time for the code increases commensurately. With care, a relative accuracy of one part in 10 4 
to 10 5 is possible. This accuracy is sufficient for confidence in the self-force values computed via matched expansions, 
presented in Sec. |VII D| 
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FIG. 17: 'Partial Field' Generated by a Static P artic le at p — 0.5. The dotted line shows the Green function. The solid line 
shows the partial field ^partial defined in Eq. (7.1). This may be interpreted as the portion of radiative field generated by 
the segment of the static-particle worldline between t' — — oo and t' — t — At. Here At = T + p* + , where T is the time 
coordinate shown on the horizontal axis. 



D. The Self-Force on a Static Particle 



In this section, we present a selection of results for a specific case: a 'static' particle at fixed spatial coordinates. 
Our goal is to compute the self-force as a function of p, to demonstrate the first practical application of the Poisson- 
Wiseman- Anderson method of matched expansions [36J [37] . 

The radiative field may be found by integrating the Green function with respect to r', where dr' = (1 — p 2 ) 1 / 2 dt f . 
Integrating a mode sum like (4.34) with respect to t' is straightforward; we simply multiply each term in the sum by 



a factor l/{iui n ). Hence it is straightforward to compute a partial field defined by 



^partial (At) = q (l - (? 



,1/2 



t-At 



G iet (t-t',p = p',i = 0)dt'. 



(7.1) 



This may be interpreted as "the field generated by the segment of the static-particle world line lying between f = —oo 
and t f = t — At" . In the limit At — > 0, the partial field ^partial will coincide with the radiative field &r. An example 
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Here, q 1 ^partial is plotted as a function of T = At 



of this calculation is shown in Fig. 
particle near spatial infinity, p — > 1. We have used the method of smoothed sums (Sec. VI D), with Z cut 



for a static 
200. The 

'partial field' ^partial shares singular points with G ret . Figure [171 shows that a significant amount of the total radiative 
field arises from the segment of the worldline after the first singularity. The Green function tends to zero in the limit 
At — » (for £ = 1/6). On the other hand, the partial field tends to a constant non-zero value in this limit. The 
constant value is the radiative field q~ 1 <&R. 

An accurate value for the total radiative field is found by using the quasilocal Green function to extend ^partial 
to coincidence, At — > + . The method is illustrated in Fig. 18 (left plot). The dashed line shows the quasilocal 
contribution to the radiative field, and the dotted line shows the distant past contribution to the radiative field, as 
a function of matching time. The former is the result of integrating from the matching point r — At to coincidence, 
and the latter from integrating from — oo to the matching point. Here, At varies linearly with the x-axis scale T (see 
caption). The right plot illustrates the same calculation for t he ra dial self- force. Here, the dotted line representing 
the contribution from the distant past is found from the sum ( |6.19 ). 

Figure [19] shows the total radiative field $>r generated by a static particle in the Nariai spacetime. The field is 
plotted as a function of p, for two cases: £ = 1/6 and £ = 1/8. In the former case, the field is negative. In the latter 
case, the field is positive, and about two orders of magnitude greater in amplitude. In both cases, the amplitude of 
the field is maximal at p = and tends to zero as p — > 1 as &r ~ (1 — p 2 ) 1 / 2 . 

The results of the matched expansion results are 



Figure [20] shows the self-force ma p acting on a static particle, 
shown as points, and the results of the 'massive field regularization method' (described in Sec. VI B) are shown as the 



solid line. The latter method provides an independent check on the accuracy of the former. We find agreement to 
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Radiative Field via Matched Expansions 



Self-Force via Matched Expansions 



'a- -0.005 



^x-*' QL 



QL+DP 



DP 



0.5 1 1.5 2 

matching time, T = t - 1' - p. - p." 




0.5 1 1.5 2 

matching time, T = t - 1' - p* - p*' 



FIG. 18: Illustration of the matching calculation of the self-force on a static particle at p — p' — 0.5. The blue dashed line 
shows the quasilocal contribution, integrated from r — At = (1 — p 2 ) 1//2 [T + p* + p*] to coincidence. The green dotted line 
shows the distant past contribution, integrated from — oo to r — At. The red solid line shows the total. In matching region 
1 < T < 2 the total approaches a constant, which corresponds to the value of the radiative field (left plot) and radial self- force 
(right plot). 




0.1 0.2 0.3 0.4 0.5 0. 




FIG. 19: The Radiative Field Generated by the Static Particle. The plot shows the radiative field for the static particle at 
p — p . For the case £ = 1/6 (left) the field is negative whereas for the case £ = 1/8 it is positive (right). Note the differing 
scales on the vertical axis. 



approximately six decimal places between the two approaches. We find the self-force at p = to be zero, as expected 
from the symmetry of the spacetime. The self- force also tends to zero as p — > 1. Between these limits, the self-force 
rises to a single peak, the magnitude and location of which depends on the curvature coupling £. We find that the 
peak of the self- force is approximately 4.9 x 10 -4 for £ = 1/6 and approximately 3.8 x 10 -2 for £ = 1/8. 



VIII. DISCUSSION AND CONCLUSIONS 



In this paper we have presented the first practical demonstration of a self-force calculation using the method of 
matched expansions. The matched expansions method was first proposed over a decade ago by Poisson and Wiseman 
[36] . We have shown that the 'quasilocal' expansion in coordinate separation [43j [46j [47j [48], valid only in the 
normal neighbourhood, may be accurately matched onto a mode sum expansion, valid in the distant past. Through 
matching the 'quasilocal' and 'distant past' expansions, the full retarded Green function may be reconstructed. With 
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FIG. 20: The Radial Self-Force on the Static Particle. The radial component of the self-force, q~ 2 ma p , is plotted as a function 
of the p of the static particle. The self-force was calculated by two methods: 'matched expansion' (black dots) and 'massive 
field regularization' (blue curve). The left plot shows the curvature coupling £ = 1/6 whereas the right plot shows £ = 1/8. 
Note the difference in scales for the two cases. 



full knowledge of the Green function, one may accurately compute the 'tail' contribution to the self-force. In this 
work, we employed the matching method to numerically compute the self- force acting on a static particle, and showed 



that the resulting self- force is in excellent agreement with the result from an alternative method (Sec. VI B), to 
approximately one part in 10 6 . 

The key new ingredient in our formulation is the so-called 'quasinormal mode sum' expansion for the distant past 
Green function. Following Leaver's approach [49], the integral over frequency in the mode sum expansion of the Green 



function (4.1 ) may be performed by deforming the contour of integration in the complex frequency plane. Poles of the 
Green function arise at (complex) quasinormal mode frequencies in the lower half-plane. The sum over the residues 
of the poles gives the quasinormal mode sum - a key contribution to the Green function (see below). 

The QNM sum is only valid at 'late times', t — t' > t c , where t c is approximately the time it takes for a geodesic 
to reflect from the peak of the potential barrier. We have demonstrated that there is a sufficient regime of overlap in 
t — t' in which both the 'quasilocal' and 'QNM sum' expansions are valid for the method to be applied successfully. 

The QNM spectra of black holes have received much attention in the last three decade s. In some studies [52, 53, 1114] . 
approximations to Q NM frequencies are found by replacing the effective potential V t with the so-called Poschl- 



Teller potential (3.7). The advantage of this replacement is that the QNMs and radial solutions of the Poschl-Teller 
potential are known in closed form. In this paper, we have taken the idea a step further. We have shown that 
the Poschl-Teller potential arises naturally if we consider the radial equation resulting from waves on an alternative 
spacetime: the product spacetime (IS2 x S 2 , first introduced by Nariai [38, 39 in 1950. This is an Einstein spacetime 
of constant scalar curvature. 

The symmetry of Nariai spacetime undoubtedly makes calculations easier. For example, geodesic motion may be 
separated into motion on the two submanifolds, (IS2 and S 2 ; the Van Vleck determinant may be written in closed form 



(5.43); and the decay rate of the quasinormal modes is independent of I. We view the Nariai spacetime as an excellent 
testing ground for our methods. Nevertheless, we should not forget the overall goal of the Poisson- Wiseman- Anderson 
proposal [36, 37 j: accurate matched-expansions calculations on physical black hole spacetimes. Below, we review some 
of the insights provided by our 'experiment' with the Nariai spacetime which will help any future calculations. 

At late times t — t' > p* + p^ it has been previously been established [67 that the quasinormal modes provide 



a complete basis on the Nariai spacetime. Here, we demonstrated that, at late times, the QNM sum ( |4.22 ) fully 
describes the retarded Green function. This is not expected to be the case on the Schwarzschild spacetime |49j. In the 
latter case, there arises a branch point at uo = and a branch cut in the frequency integral which gives a 'power-law 
tail' contribution [50] to the Green function (Fig. [6j. The branch point is notably absent for the Nariai spacetime, 
making the analysis simpler. In this paper we have thoroughly investigated the effect of the quasinormal modes. The 
contribution of the branch cut integral to the self-force remains to be quantified. We hope to pursue this calculation 
in a forthcoming work. 

This study has provided a number of insights into the properties of the Green function, which are of relevance to 
any future investigation of the Schwarzschild spacetime. Namely, 

• The Green function G re t(x, x') is singular whenever x and x' are connected by a null geodesic. The nature of the 
singularity depends on the number of caustics that the wave front has passed through. After an even number 
of caustics, the singularity is a delta-distribution, with support only on the light cone. After an odd number 
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of caustics, the Green function diverges as l/Va, where a is the Synge world function. A four- fold repeating 
pattern occurs, i.e. 5, l/ira, —5, —l/ira, 5, etc. 



The four-fold singular structure can be shown to arise from a Hadamard-like ansatz (5.38) valid even outside 
the normal neighborhood, if we allow U(x,x') to pick up a phase of —i upon passing through a caustic. The 
accumulation of phase may be deduced by analytically continuing the integral for the Van Vleck determinant 
through the singularities (due to caustics). 



Hadamard's form for the Green function (4.35) is only strictly valid if x and x' are in a convex normal neighbour- 
hood [41]. The extension beyond the normal neighbourhood does not seem to be known. However, this work 
and some previous studies jl 15] would tentatively suggest that, if x and x' are connected by a countable number 
of distinct geodesies, then the Green function may be found from the sum of their Hadamard contributions. 



• The effect of caustics on wave propagation has been well-studied in a number of other fields, such as optics |116| , 
acoustics [93] , seismology [83 , symplectic geometry [94] and quantum mechanics [95]. It may be that mathe- 
matical results developed in other fields may be usefully applied to wave propagation in gravitational physics. 

A further observation made in this work is confirmation that the 'tail' self-force cannot be calculated from the 



'quasilocal' contribution alone. For instance, Fig. [TT] would appear to show that a significant part of the radiative 
field is generated by the segment of the world line which lies outside of the normal neighbourhood (i.e. beyond the 
first caustic). Unlike in flat space, the radiated field generated by an accelerated particle may propagate once, twice, 
etc. around the black hole before later re-intersecting the world-line of the particle (see Fig. [2|. Radiation from near 
these orbits will give an important contribution to the self-force which cannot be neglected. 

Let us conclude by examining the prospects for a practical 'matched expansions' calculation on the Schwarzschild 
spacetime. Happily, the 'quasilocal' expansion is now in excellent shape, as described in [47 ; the challenge remains 
the 'distant past' expansion. We have already mentioned that a quasinormal mode sum will not be sufficient; it must 
be augmented by a branch cut integral. We have reasons to be optimistic that this is a tractable calculation [50] . 
Perhaps more difficult will be the accurate numerical computation of QNM frequencies and radial functions. A further 
difficulty will be in integrating the mode sum over the worldline; for accuracy we wish to avoid numerical integration if 



possible. For the static particle, it was straightforward to integrate each term in the mode sum analytically (Sec. VI C ) 
with respect to time. For other trajectories (e.g. circular orbits) this may present more of a challenge. Finally, we note 
that a range of established results are available on Schwarzschild. For the static particle the self-force is zero [100 , 117 ; 
for radial trajectories and circular and eccentric orbits, accurate numerical results are available [26] [29] I102[ I118j . 
This will surely help the validation of the 'matched expansions' method. We hope to undertake such a study in the 
near future. 
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APPENDIX A: LARGE-A ASYMPTOTICS OF 2 Fi(§ - zA, | — iX; 1 - 2i\; -e~ T ) 

To determine the singularity structure of the Green function in Sec. |V A| we required the large A asymptotic behaviour 
of 2^1 (— /5, — /5; —2/3; z) where (3 = — \ + iA, z — —e~ T . The required asymptotics may be found by applying the 
WKB method jl 19] to the hypergeometric differential equation 

z(l - - [2/3(1 -z) + z]f z -(3 2 u = (Al) 

which has solutions u(z) = 2 Fi(-(3, -/3; -2/3; z). Inserting the WKB ansatz u(z) ~ e *So(z)+s 1 (z)+\- 1 s 2 (z)... j mmec ii- 
ately yields a quadratic equation for Sq, 



z(l - z){S'o) 2 - 2i(l - z)S' + 1 = 0. 



(A2) 
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In our case, z = — e T , we require the root which is finite at z = 0. We impose 5q(0) = to get 



In 



-T - 2 In 2 + In 



s'^mi-z))- 1 



Si 



ln(l + e- T ) . 



4 [yT^-ljy ~ " V nXi f - 

At next order, we obtain the equation 

2 [i(l - z) - z(l - z)^] ^ = z(l - *)Sff + (1 - 2z)S' + i. 
It is straightforward to show that this reduces to 

1 



(A3) 



(A4) 



(A5) 



Inserting (A5) and (A3) into ffl) leads to the quoted result, Eq. (5.2). Of course, the asymptotic approximation may 



be further refined by taking the WKB method to higher orders. 



APPENDIX B: POISSON SUM ASYMPTOTICS 

In this appendix w e deri ve asymptotic approximations for t he sin gular structure of the Green function using the 
Poisson sum formula (5.18). Our starting point is expression (5.29) for the 'n=0' fundamental modes, in which the 



Legendre polynomials P^(cos7) have been replaced by angular waves Q^ 1 ^ 2 ( cos l) which are, in turn, approximated 



by Hankel functions Hq^\i/j) using (5.28). 



Let us consider the X\ and Z 2 integrals arising from substituting (|5.29 [) into (|5.18|). These integrals are singular at 

(Bl) 



X = 27T — 7 and \ = 2tt + 7, respectively. First, let us consider X\ ( |5.18[ ) which can be written 



.4( 7 )Re / dz/(-w)i/ 2 e i(x - 27r ^^ +) (i/ 7 ) 



with .4.(7) as defined in Eq. (5.32) . 

For x > 27T — 7, the integral may be computed by rotating the contour onto the positive imaginary axis (y = iz) to 
obtain 



2.4(7) 



dzz 1/2 e- (x - 2 ^ z K (jz) 



^(3/2, 1/2; 2; 



2[x-(2tt-7)] 3 / 2 



(B2) 



Here we have applied the identity Hq + \zx) — 2Kq(x) / '(in) , where Kq is the modified Bessel function of the second 
kind, and the integral is found from Eq. 6.621(3) of Ref. [65] . 

For x < 27r — 7, the integral may be computed by rotating the contour onto the negative imaginary axis (v = — iz). 
First, we make the replacement (uj) = 2J (vj) — Hq~\vj) and note H^~\—ix) = 2K (x)/(— in) to obtain 



Ti « ^^Re J dzz 1 ^e-^-^ z [7rI (jz)^iK (jz)} 



(B3) 



Here Iq is a modified Bessel function of the first kind. Since we are taking the real part, the Kq term is eliminated, 
and we obtain 



2.4(7) 



(2tt-7-x) 1 (27t + 7-x)' 



-1/2 



E 



27 



2tt + 7 - x 



(B4) 



where E is the elliptic integral of the second kind defined in, for example, Eq. 8.111(3) of Ref. [65] . 

The X2 integral may be calculated in a similar manner. For x < 27r + 7, we rotate the contour onto the negative 
imaginary axis, 



i(x-27i> 



I2 « -^(7)Re / dvi-ivf^H^iv^y 
Jo 

~ Mkll Re i f°° dzz^e-^-^Ko^z) = 

7T Jo 



(B5) 
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For x > 27r + 7, we rotate the contour onto the positive imaginary axis after taking the complex conjugate 

POO 

Jo 

~ ^4(7) Re P ^V2 e -(x-^)^ (i7r/ ( 7 ^) + #0(7*)) (B6) 



o 



The imaginary term does not contribute and hence X2 is equal and opposite to T\ defined by Eq. (B2 ) when \ > 27T+7. 



APPENDIX C: GREEN FUNCTION ON T x § 

To the best of our knowledge, the four-fold singularity structure for the Green function of Sec. [V] has not been 
shown before in the literature (with the exception of [84j) within the theory of General Relativity. We therefore wish 
to illustrate its derivation and manifestation in the simplest of spacetimes including S 2 -topology: 

ds 2 = -dt 2 +dQl, (CI) 

where Synge's world function is simply given by a = \ [—At 2 + A(/> 2 ), in the case of a conformally-coupled (£ = 1/8) 
scalar field. 

Let x = (t, 0, (j)) denote any point in this spacetime. We introduce the Wightman function G+(x,x') (it satisfies 
the homogeneous scalar wave equation - see, e.g., [86]), from which the 'retarded' Green function G ret (x,x f ) is easily 
obtained: 

+00 +Z A +00 



G + (x, x ') = j2 £ ^w*L(^) = ^E e " i( ' +1/2)Atp '( cos ^)' ( C2 ) 

1=0 m=-l 1=0 

G ret {x,x') = -20(At) lm(G + (x,x')), (C3) 



where At = t — t' , <&i m (x) = e ^ +1 / 2 ^l^ m (#, 4>)/ \J (21 + 1) are the Fourier-decomposed scalar field modes on (CI) 
normalized with respect to the scalar product 

(<S>im^i>m>) = ~i J dVn» [^ lrn (x)d^ rnf (x) - ^l rnf (x)d^i rn (x)} = 5 w 5 mm ,, (C4) 

where £ is a Cauchy hypersurface with future-directed unit normal vector n M and volume element dV. 



We now apply exactly the same tricks as in section VP in order to derive the four-fold singularity structure in the 



Green function from the large-/ asymptotics of the field modes. We use the Poisson sum formula 

+00 +00 r+00 



£>(Z + l/2)= J2 (-W dvg{v)e 2 ™» (C5) 

1=0 s=-oo '' {) 



to re- write the mode sum in (C2) as 

4nG + (x,x') = £ ("I)* / ^e-- A *P [/ _ 1/2 (cos 7 )e 2 ' ri - = £ G», 

s=-oo ^° N=l 



(C6) 



r + oo 

G%{x,x')= / ^Ar(cos 7 )e- ii/A£ . (C7) 
Jo 

The Legendre functions P /i (cos 7 ) and Q /i (cos 7 ), as well as i?Ar(cos 7 ), are standing waves. This is in contrast to 

Q^(cos 7 ), which are travelling waves. 

We can now use large-order uniform asymptotics (see [85 ], I120| ) for the Legendre functions: 

x 1/2 

7 x 



P i/ _ 1 / 2 (cos 7 ) ~ ^— J Jo(z/ 7 ), — > 00, "valid in a closed uniform interval containing 7 = 0", (C8) 

/ \ 1/2 

Qi/-i/2(cos7) - ( J IoC^t), M -> "valid with respect to 7 e (0,tt/2]", (C9) 

^ \ sin 7 y 

QS /2 (cos 7 ) - \ (£^) 1/2 ^o T) (^), M - 00, (CIO) 
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to leading order. 

The contribution to the Wightman function from the TV = 1 orbit wave is immediately obtained by using the large- 
order asymptotics of the Legendre function P^(cos7) only, which are "valid in a closed uniform interval containing 
7 = 0"- this is wh at w e will mean by a result being valid "near" 7 = 0. Similarly, we can obtain a result valid "near" 
7 = 7r by using in (C2) the symmetry P^(cos7) = (— l)'Pf (cos(7r — 7)) for I G N. We then obtain for N = 1: 



AttG^ 1 (x,x') 



7 



1 



sin 7 ^72 _ At 2 



7T — 7 



7T — 7 



(tt - 7) V(^-7) 2 -(At-^) 2 V sin(^ - 7) y/-(At-*y) [At - (2tt - 7)] ' 



"near" 7 = 
(Cll) 



near 7 = 7r 



(C12) 



where, for convergence, a small imaginary part was given to At and/or 7, in agreement with the Feynman prescription 
L (j — > (j -\- ie\ The result for G^ =1 (x, x') valid "near" 7 = is singular at At = ±7, corresponding to a = b efore a 
caustic has been crossed. It is in accord with the Hadamard form in 3-D [121 and the Van Vleck determinant (5.48), 
before a caustic has been crossed (and so without the phase factor). The result for G+ =1 (x,x f ) valid "near" 7 = tt 
is singular at At = 7, corresponding to the case where it has not gone through any caustics, and at At = 2tt — 7, 
corresponding to the case where it has gone through one caustic; it has thus picked up a factor "— i", as expected. 
Note that these zeros inside the squared root in the denominator are simple zeros along the null geodesic, except at 
the caustic point itself, where the two zeros coincide and so i t bec omes a double zero. 
Similarly to N = 1, we can use (C6) and the asymptotics (CIO) together with [65] 



I±{T,l) 



1 



Vt 2 - (t - W 



1 T ^ In (iX + v 7 ! -X 2 ) 



T e R, e > |Im 7 |, (C13) 



where X = (T — ie)/j (again, a small imaginary part needs to be given for convergence, in accordance with the 
Feynman prescription), in order to obtain for N > 1: 



5 {-l) N ' 2 J [7+ (At + Nit, S) + J_ (At — Nir, 5)] , "near" 7 : 

sinS 2 1 -i [7+(Ai + (AT - 1)tt, H) + 7_ (Ai - (N + 1)tt, H)] , "near" 7 : 







(C14) 



for N even, where 5 = 7 "near" 7 = and S = 7r — 7 "near" 7 = 7r. For T V odd, merely: (1) swap I± — > 7 T , and (2) 
replace by N - 1 if "near" 7 = or N by iY + 1 if "near" 7 = tt in (|C14|). We can re- write: 7 2 - (At ± Nn - ie) 2 = 



• [(At - ie) - (tNtt - 7)] [(At - ie) - (7 T Nit)]. Note, however, that I±(T, 7) is regular at X = =pl. 
We then have that the singular behaviour goes as 



AttG + (x,x') 



r+/+(At, 7 ), 



< At < tt, 7-O 

1 / S J — U+(At — 7T, 7T — 7), 7T < At < 27T, 7 ~ 7T 



2 V sinS j -i + (At-27r,7), 2tt < At < 3tt, 7 - 

[+ii+(At - 37T,7T - 7), 37T < At < 47T, 7 - 7T 



(C15) 



The 4- fold singularity structure arises clearly: a phase of tt/2 is picked up everytime the null geodesic joining x and 
x' goes through a caustic (7 = or tt). 

The expression for G + (x,x / ) is simplified by noting that, "near" 7 = 0: 



4tt [G%(x,x') + G% +1 (x,x')] 
1 



7 
sin 7 



V7 2 - (At — Ntt — ie) 2 V7 2 - (At + Ntt - ie) 2 



N even 



7T — 7 

sin(7r — 7) 



^(tt - 7 ) 2 - (At - Ntt - ie) 2 y/(7r - 7 ) 2 - (At + Ntt - ief 



N odd. 



(C16) 
(C17) 

(C18) 



Similarly "near" 7 = 7r. The Poisson sum formula has yielded a sum over geodesic paths, labelled by the index iV, 
and the large-order asymptotics for the Legendre functions have yielded the correct singularity structure near the null 
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geodesies, allowing for the correct phase change at each caustic. 
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